//
//  ASolverUnitTest.h
//  pbsolvers
//
//  Created by Marielle Soniat on 10/6/15.
//  Copyright (c) 2015 Marielle Soniat. All rights reserved.
//

#ifndef pbsolvers_ASolverUnitTest_h
#define pbsolvers_ASolverUnitTest_h

#include "ASolver.h"

class ASolverUTest : public ::testing::Test
{
public :
  
protected :
  
  int vals_;
  shared_ptr<Constants> const_;
  vector< shared_ptr<BaseMolecule> > mol_;
  vector< shared_ptr<BaseMolecule> > mol_sing_;
  
  virtual void SetUp()
  {
    const_ = make_shared<Constants>();
    
    mol_.clear( );
    Pt pos[2]     = { Pt( 0.0, 0.0, -5.0 ), Pt( 10.0, 7.8, 25.0 ) };
    Pt cgPos[2]   = { Pt( 0.0, 0.0, -5.5 ), Pt( 11.0, 6.9, 24.3 ) };
    double cg[2] = { 5.0, -0.4}; double rd[2] = { 5.6, 10.4};
    
    Pt cgPosSi[2] = { Pt( 0.0, 0.0, -35.0 ), Pt( 0.0, 0.0, 0.0 ) };
    
    shared_ptr<MoleculeAM> molNew, molSing;
    
    for (int molInd = 0; molInd < 2; molInd ++ )
    {
      int M = 1;
      vector<double> charges(M);
      vector<double> vdW(M); vector<Pt> posCharges(M);
      charges[0] = cg[molInd]; posCharges[0] = cgPos[molInd]; vdW[0] = 0.0;
      
      molNew = make_shared<MoleculeAM>("stat",rd[molInd],charges,posCharges,
                                       vdW,pos[molInd],
                                       molInd, 0);
      mol_.push_back( molNew );
      
      charges[0]    = 2.0; posCharges[0] = cgPosSi[molInd];

      molSing = make_shared<MoleculeAM>( "stat", 10.0, charges, posCharges, vdW, molInd, 0);
      mol_sing_.push_back( molSing );
    }
  } // end SetUp
  
  virtual void TearDown() {}
  
  double ATrip0[15]   = {6.01094069,0.0141280225,1.03433512,-0.804854118,
    1.24302352e-05,2.58107532e-05,5.8366607e-05,-0.182460812,3.82002752e-06,
    0.235556808,0.162128207,-7.8270901e-09,3.91410002e-07,
    2.55495441e-08,0.226073754};
  double ATrip0im[15] = {0,0,1.03698189,0,0.000412418784,3.23186115e-06,0,
    -0.18242449,3.11917432e-07,-0.235553474,0,2.63554468e-06,1.04180458e-08,
    5.16522503e-07,1.02665424e-08};
  double ATrip1[15]   = {6.01119603,-0.00735850845,1.01433618,-0.80724428,
    0.000252370073,0.00231355778,-3.33699295e-05,-0.182293154,-2.53839549e-06,
    0.23531127,0.16214575,1.7126064e-06,-1.9023086e-05,-9.1676707e-08,
    0.22609914};
  double ATrip1im[15] = {0,0,1.03337432,0,0.000193710438,-4.13320113e-05,0,
    -0.182483415,-9.28886031e-06,-0.235542082,0,1.32479048e-06,1.79437885e-06,
    2.38805876e-07,-2.11951326e-06};
  double ATrip2[15]   = {6.01116143,-0.00712987242,1.05234476,-0.807013253,
    -0.000264824695,0.0019962685,-3.00198469e-05,-0.182594503,-3.60529882e-06,
    0.235757586,0.162141472,-1.7044803e-06,-1.46864088e-05,6.58681972e-08,
    0.22609354};
  double ATrip2im[15] = {0,0,1.03341635,0,0.000184134752,3.80938095e-05,0,
    -0.182482137,9.60094345e-06,-0.235542368,0,1.19249837e-06,-1.80475738e-06,
    2.8073034e-07,2.10915172e-06};
  
  double ATripPBC3[15] = {6.01060565,1.7911869e-16,1.03422506,-0.805960634,6.04422072e-17,0,-7.74056204e-17,-0.182461408,0,0.235556664,0.162125951,-2.21982325e-17,0,1.95770009e-17,0.226073838};
  double ATripPBC3im[15] = {0,0,1.03422506,0,6.04422072e-17,6.04422072e-17,0,-0.182461408,4.32675923e-33,-0.235556664,0,-2.21982325e-17,-1.04643472e-17,-1.95770009e-17,-2.76860602e-17};
  
  double APBCTrip0[15]   = {6.01105542,-0.00074451012,1.03428935,-0.80430679,0.000926503077,7.23059366e-05,-1.26250549e-05,-0.182464792,-2.83937838e-07,0.235556133,0.162128461,6.99600862e-06,1.36983475e-06,-7.4788967e-07,0.226073548};
  double APBCTrip0im[15] = {};
  double APBCTrip1[15]   = {6.01083506,-0.0103527557,1.03144083,-0.805073567,0.000467065091,3.623056e-05,-4.96995897e-05,-0.182508586,-6.13429469e-06,0.23555845,0.162127705,3.76497861e-06,7.18712557e-07,-3.3326665e-07,0.226073712};
  double APBCTrip1im[15] = {};
  double APBCTrip2[15]   = {6.01082591,0.00960836601,1.03707364,-0.80519449,0.000459081726,3.60474528e-05,3.70745736e-05,-0.182417611,5.8505518e-06,0.235554346,0.162126705,3.22663862e-06,6.50256324e-07,-4.14147247e-07,0.226073674};
  double APBCTrip2im[15] = {};
  
  double APBCoutside[15]   = {6.01060565,1.7911869e-16,1.03422506,-0.805960634,
    6.04422072e-17,0,-7.74056204e-17,-0.182461408,0,0.235556664,0.162125951,
    -2.21982325e-17,0,1.95770009e-17,0.226073838};
  double APBCoutsideim[15] = {0,0,1.03422506,0,6.04422072e-17,6.04422072e-17,0,
    -0.182461408,0,-0.235556664,0,-2.21982325e-17,-1.04643472e-17,
    -1.95770009e-17,-2.76860602e-17};

  double A0[15] = {5.06502332,-0.459077988,-0.000488586514,0.0313638895,
    -0.0001113252,-2.2608511e-05,-0.00206626713,-1.20922746e-05,-5.8562079e-06,
    -9.08691077e-07,0.000132860057,-3.86433215e-07,-9.97731223e-07,
    -3.82176443e-07,-3.75510255e-08};
  double A0_im[15] = {0,0,-0.00021783804,0,-3.85223997e-05,-2.08028337e-05,0,
    -8.64506194e-07,-2.89571751e-06,-1.2658108e-06,0,9.07199987e-07,
    2.00374768e-07,-2.07633434e-07,-8.16687899e-08};
  double A1[15] = {-0.396665696,-0.157791465,-0.100108147,0.113721012,
    0.0563843732,0.00185962096,-0.0415378387,-0.0286835422,-0.00276725208,
    0.00189910059,0.0110646934,0.013544291,0.00205846384,-0.00116220012,
    -5.31391954e-06};
  double A1_im[15] = {0,0,0.00922544874,0,0.0293832761,0.0232705597,0,
    -0.0222273022,-0.0142638165,-0.000444043694,0,0.0109506821,0.00839150215,
    0.00171418682,0.000215801105};
  
  double AMul0Sing[15]   = {6.0131921,0.451572318,1.06409876,-0.565532819,
    0.0265785191,0,0.107782872,-0.166377441,0,0.235821288,0.208159863,
    0.00831875258,0,0.000298201925,0.226117001};
  double AMul0Singim[15] = {0,0,1.06409876,0,0.0265785191,0,0,-0.166377441,0,
    -0.235821288,0,0.00831875258,0,-0.000298201925,0};
  double AMul1Sing[15]   = {6.0131921,-0.451572318,1.06409876,-0.565532819,
    -0.0265785191,0,-0.107782872,-0.166377441,0,0.235821288,0.208159863,
    -0.00831875258,0,-0.000298201925,0.226117001};
  double AMul1Singim[15] = {0,0,1.06409876,0,-0.0265785191,0,0,-0.166377441,0,
    -0.235821288,0,-0.00831875258,0,0.000298201925,0};
  
  double AMul0SingF[15]   = {6.01581722,0,1.09432855,-0.315346337,0,0,0,-0.150033773,0,0.236085919,0.257440018,0,0,0,0.226160163};
  double AMul0SingimF[15] = {0,0,1.09432855,0,0,0,0,-0.150033773,0,-0.236085919,0,0,0,0,0};
  double AMul1SingF[15]   = {6.01423317,0.542415383,1.06836766,-0.549177154,0.028691067,0,0.107338231,-0.165655581,0,0.235823283,0.206234945,0.00853035235,0,0.00029923875,0.226117093};
  double AMul1SingimF[15] = {0,0,1.06836766,0,0.028691067,0,0,-0.165655581,0,-0.235823283,0,0.00853035235,0,-0.00029923875,0};
  double AMul2SingF[15]   = {6.01423317,-0.542415383,1.06836766,-0.549177154,-0.028691067,0,-0.107338231,-0.165655581,0,0.235823283,0.206234945,-0.00853035235,0,-0.00029923875,0.226117093};
  double AMul2SingimF[15] = {0,0,1.06836766,0,-0.028691067,0,0,-0.165655581,0,-0.235823283,0,-0.00853035235,0,0.00029923875,0};
  
  double A0Sing[15] = {2.08493624,0.0595104925,0,0.0264609804,0,0,
    0.00906707868,0,0,0,0.00287107389,0,0,0,0};
  double A1Sing[15] = {2.08493623, -0.059510451, 0, 0.026460955, 0, 0,
    -0.00906706693, 0, 0, 0, 0.00287106903, 0, 0, 0, 0};
  
  double dTA00X[15]   = {8.40366148e-05,1.24371239e-05,-0.00024223864,
    7.32052642e-07,-2.06649326e-05,1.55421909e-08,-1.33035844e-08,
    -8.90141524e-07,3.28549285e-08,-1.99245602e-07,-8.66506018e-09,
    8.76895226e-09,6.77273003e-09,-3.3842657e-08,-1.12365334e-09};
  double dTA00Xim[15] = {0,0,2.28653741e-06,0,2.54344891e-07,-5.70164427e-06,0,
    9.68136519e-09,-4.93159502e-07,1.10213427e-08,0,-1.48696517e-09,
    -1.63296663e-08,3.05594515e-09,-1.47030724e-08};
  double dTA00Y[15]   = {0.00210486595,0.000412647342,2.28653741e-06,
    4.61199423e-05,2.54344891e-07,1.2074341e-05,3.99665825e-06,9.68136519e-09,
    1.64053954e-06,1.53123411e-08,2.88252398e-07,-1.48696517e-09,
    1.63396968e-07,1.89505581e-09,1.2600396e-08};
  double dTA00Yim[15] = {0,0,-0.000278760786,0,-2.68948152e-05,-2.72384346e-07,
    0,-1.62340475e-06,-2.38996886e-08,1.55756718e-08,0,-6.15154603e-08,
    -1.25082386e-09,7.72893677e-09,-3.70533308e-10};
  double dTA00Z[15]   = {0.00762745576,0.00098754963,8.79437465e-06,
    6.75899726e-05,6.39735563e-07,3.11494127e-06,3.23452259e-06,-8.88300647e-09,
    4.6722108e-07,2.31699754e-08,7.01440824e-08,-7.6619317e-09,
    5.01207285e-08,4.29919861e-09,-1.46977779e-08};
  double dTA00Zim[15] = {0,0,0.000291785734,0,3.87091775e-05,2.54344891e-07,0,
    3.54874075e-06,1.26798933e-08,4.68415939e-07,0,2.63294193e-07,
    -2.06841692e-09,7.88415973e-08,1.75044317e-09};
  
  double dTA11X[15]   = {-0.0151853099,0.000252509934,0.00211174413,
    0.000216265654,1.16422283e-05,-0.000293015978,2.59871423e-06,
    -2.97749845e-05,-7.96572085e-07,3.73320959e-05,-3.46401015e-06,
    6.28382238e-09,3.50804631e-06,1.54270798e-08,-4.49835015e-06};
  double dTA11Xim[15] = {0,0,-2.92423424e-05,0,-7.57435756e-06,2.00432626e-05,
    0,1.68562581e-06,-2.72333015e-07,-3.1580328e-06,0,-9.15146217e-08,
    -1.55301981e-07,2.94339673e-08,4.47759661e-07};
  double dTA11Y[15]   = {-0.000649545295,0.000193817791,-2.92423424e-05,
    -2.76416645e-05,-7.57435756e-06,9.08336352e-06,2.0103999e-06,1.68562581e-06,
    8.0280573e-07,-2.27851923e-06,-8.59235284e-08,-9.15146217e-08,
    -2.08106667e-07,-3.86472161e-08,3.7684397e-07};
  double dTA11Yim[15] = {0,0,-0.00116193389,0,1.57819644e-05,0.00019713346,0,
    6.07902481e-06,-5.9292554e-07,-2.75674893e-05,0,5.2794759e-08,
    -1.30545444e-06,1.85252068e-09,3.55102393e-06};
  double dTA11Z[15]   = {-0.00397272143,-0.000901392458,0.000178551487,
    -3.84530678e-05,0.000168987733,-2.06986804e-06,2.69764415e-05,
    2.30411642e-06,-2.23655072e-05,-8.31383539e-08,-7.24390524e-08,
    -3.07335924e-06,-3.31687119e-08,2.55517124e-06,4.79933141e-09};
  double dTA11Zim[15] = {0,0,0.000137049874,0,-2.28940877e-05,-7.57435756e-06,
    0,1.78217115e-06,2.12661856e-06,2.16564579e-07,0,-8.01151965e-08,
    -1.30447713e-07,-1.93195131e-07,-3.25737536e-09};
  
  double dTA22X[15]   = {0.0138345289,-0.000264971459,0.00181443092,
    -0.000171521013,9.31304576e-06,0.000239127813,-2.5849498e-06,
    -2.31391192e-05,7.63301685e-07,2.90615005e-05,2.61981372e-06,
    -1.4916962e-08,-2.62055104e-06,1.90680347e-08,3.33215071e-06};
  double dTA22Xim[15] = {0,0,2.69513189e-05,0,7.82883757e-06,
    1.94376751e-05,0,-1.6952699e-06,-1.71150792e-07,3.1468793e-06,0,
    9.0024276e-08,-1.67087182e-07,-2.63740263e-08,4.4839425e-07};
  double dTA22Y[15]   = {-0.000617451849,0.000184236798,2.69513189e-05,
    -2.60347135e-05,7.82883757e-06,9.11776478e-06,1.80861763e-06,
    -1.6952699e-06,7.94752641e-07,2.26309537e-06,-6.50883031e-08,
    9.0024276e-08,-2.06803856e-07,4.05436927e-08,3.7554146e-07};
  double dTA22Yim[15] = {0,0,-0.00101028362,0,1.51927385e-05,-0.00016218759,0,
    4.55402329e-06,6.16985011e-07,-2.14735717e-05,0,4.94370144e-08,
    9.54333145e-07,4.82955894e-09,-2.62572648e-06};
  double dTA22Z[15]   = {-0.0038492851,-0.000721285509,-0.000187363115,
    -3.47275886e-05,-0.000133129258,-2.93984638e-06,2.11442689e-05,
    -2.29483949e-06,-1.7246713e-05,5.97335319e-08,-4.43405578e-08,
    2.32349965e-06,-4.59035476e-08,-1.89632039e-06,5.03406139e-09};
  double dTA22Zim[15] = {0,0,0.000130275089,0,-2.15659715e-05,7.82883757e-06,0,
    1.6053902e-06,-2.13925215e-06,2.54584389e-07,0,-6.1290081e-08,
    1.28374534e-07,-1.98818893e-07,5.00973358e-09};
  
  double dTA02X[15]   = {-0.00134312154,0.000264445158,0.000116010163,
    -2.97402876e-05,-9.40166304e-06,6.18828082e-06,2.59867123e-06,
    3.26725526e-07,-7.75828158e-07,1.28879662e-07,-1.89570995e-07,
    1.67712122e-08,6.96531629e-08,-2.07034721e-08,-3.0608302e-10};
  double dTA02Xim[15] = {0,0,4.55809629e-05,0,-7.82860105e-06,-2.27082704e-06,
    0,9.29391901e-07,1.71215149e-07,-1.96522678e-07,0,-9.00289095e-08,
    -9.13508283e-10,2.63862862e-08,-7.31276472e-09};
  double dTA02Y[15]   = {0.000936010711,-0.000184235667,4.55809629e-05,
    2.07108202e-05,-7.82860105e-06,5.71633817e-06,-1.80864707e-06,
    9.29391901e-07,-7.94704561e-07,2.69290902e-07,1.31834659e-07,
    -9.00289095e-08,8.13075048e-08,-4.05349697e-08,7.61554053e-09};
  double dTA02Yim[15] = {0,0,0.00014957279,0,-1.51711926e-05,5.28159952e-06,
    0,1.01243297e-06,-6.12044688e-07,4.296989e-08,0,-4.97434464e-08,
    4.8578427e-08,-4.01876294e-09,-3.60700385e-09};
  double dTA02Z[15]   = {0.00384765899,-0.000501599176,0.000186990964,
    3.48077424e-05,-2.49543843e-05,2.88476479e-06,-1.71518864e-06,
    2.30661096e-06,-4.36863404e-07,-6.68643219e-08,4.26242074e-08,
    -1.73086554e-07,4.7423916e-08,1.129776e-08,-5.89893548e-09};
  double dTA02Zim[15] = {0,0,-0.000130274289,0,1.73783056e-05,
    -7.82860105e-06,0,-1.60541546e-06,1.18424623e-06,-2.54538487e-07,
    0,1.20374152e-07,-1.28381054e-07,4.30971946e-08,-5.002315e-09};
  
  double dTA10X[15]   = {-0.00130993449,-0.000254679919,0.000113852059,
    -2.81182583e-05,9.17038656e-06,6.22042771e-06,-2.39629293e-06,
    3.19878772e-07,7.83716188e-07,1.13261912e-07,-1.68775884e-07,
    -1.50020177e-08,7.09628726e-08,1.87642858e-08,-1.60443988e-09};
  double dTA10Xim[15] = {0,0,-4.96111547e-05,0,-8.41478095e-06,2.5486814e-06,
    0,-9.83275685e-07,1.96040513e-07,2.18428122e-07,0,-9.33933949e-08,
    2.27461873e-09,2.93631643e-08,7.6774407e-09};
  double dTA10Y[15]   = {-0.00103640083,-0.000201539018,-4.96111547e-05,
    -2.22575405e-05,-8.41478095e-06,-6.01973614e-06,-1.89755251e-06,
    -9.83275685e-07,-8.16253293e-07,-2.73962252e-07,-1.33718702e-07,
    -9.33933949e-08,-8.10510319e-08,-4.00164781e-08,-6.94320742e-09};
  double dTA10Yim[15] = {0,0,0.000137353498,0,1.3152869e-05,4.60385499e-06,
    0,7.84698074e-07,4.94619605e-07,3.11319409e-09,0,2.90853856e-08,
    3.4206926e-08,-2.38930338e-09,-5.21757356e-09};
  double dTA10Z[15]   = {-0.00375206577,-0.000479482777,-0.000180085897,
    -3.19354823e-05,-2.36134099e-05,-1.99124122e-06,-1.43640969e-06,
    -2.12927475e-06,-2.96207016e-07,1.18023186e-07,-2.101148e-08,
    -1.54295517e-07,-3.14443593e-08,1.97054319e-08,7.47892317e-09};
  double dTA10Zim[15] = {0,0,-0.000142509607,0,-1.86913911e-05,
    -8.41478095e-06,0,-1.68608032e-06,-1.25316907e-06,-2.53200807e-07,0,
    -1.22243709e-07,-1.33218596e-07,-4.22334851e-08,-3.76651522e-09};
  
  double dTA21X[15]   = {0.0139262899,-5.10438415e-07,-0.00224230427,
    -0.000245273608,8.5965464e-08,0.000299219901,1.33139083e-08,2.91976279e-05,
    -1.21544176e-08,-3.7433673e-05,3.26576252e-06,-1.79959726e-09,
    -3.4316127e-06,1.58722952e-09,4.49691933e-06};
  double dTA21Xim[15] = {0,0,7.25351982e-05,0,2.25090129e-10,-1.71675376e-05,
    0,-7.65908567e-07,-6.12575553e-11,2.9504751e-06,0,-4.41112622e-12,
    1.68007384e-07,1.16715163e-11,-4.41099212e-07};
  double dTA21Y[15]   = {-0.000318571689,-1.07581356e-09,7.25351982e-05,
    5.3241063e-06,2.25090129e-10,-1.48346991e-05,2.80247377e-11,
    -7.65908567e-07,-4.5763443e-11,2.53248806e-06,-6.67490054e-08,
    -4.41112622e-12,1.25501358e-07,8.30440738e-12,-3.83172417e-07};
  double dTA21Yim[15] = {0,0,0.00100830138,0,-2.09013662e-08,-0.000192124427,
    0,-7.14553469e-06,4.79360922e-09,2.75463418e-05,0,2.97400437e-10,
    1.35278166e-06,-7.86897492e-10,-3.55500189e-06};
  double dTA21Z[15]   = {1.57677268e-06,0.00143160944,-3.60934464e-07,
    -7.7755832e-08,-0.00019330212,5.34334151e-08,-2.48991735e-05,
    1.14219056e-08,2.26770637e-05,-6.91900298e-09,1.66576792e-09,
    2.89260044e-06,-1.47555988e-09,-2.53956994e-06,8.39380656e-10};
  double dTA21Zim[15] = {0,0,-7.60715061e-10,0,4.18783338e-06,2.25090129e-10,
    0,2.40419269e-11,-9.55044034e-07,-4.36911396e-11,0,-5.90864162e-08,
    -6.20776981e-12,1.55727906e-07,7.06255555e-12};
  
  double dTA00XSin[15]   = {0.0458894702,0,-0.140787317,0.0409529396,0,-0.00781945221,0,-0.0565498267,0,0.000110712135,0.0143359569,0,-0.00392004348,0,-0.000115663939};
  double dTA00XimSin[15] = {0,0,0,0,0,-0.00887657132,0,0,0,0,0,0,-0.00514454534,0,0.000115663939};
  double dTA00YSin[15]   = {0.0458894702,0,0,0.0409529396,0,0.00781945221,0,0,0,0,0.0143359569,0,0.00392004348,0,-0.000115663939};
  double dTA00YimSin[15] = {0,0,-0.140787317,0,0,-0.00887657132,0,-0.0565498267,0,-0.000110712135,0,0,-0.00514454534,0,-0.000115663939};
  double dTA00ZSin[15]   = {0,0.402107013,0,0,0.0334608793,0,0.130743492,0,0,0,0,0.0128272774,0,0.000654562068,0};
  double dTA00ZimSin[15] = {0,0,0,0,0.0334608793,0,0,0,0,0,0,0.0128272774,0,-0.000654562068,0};
  
  double dTA11XSin[15]   = {0.0260681318,0.0287069673,-0.0735281909,0.02122356,-0.0478521949,-0.00405945595,0.0130730001,-0.0261676531,-0.00324590507,5.5474199e-05,0.00718144537,-0.0134519833,-0.00196385203,7.54365475e-05,-5.79088536e-05};
  double dTA11XimSin[15] = {0,0,0,0,0,-0.00459199244,0,0,-0.00391062067,0,0,0,-0.00257691824,0,5.79088536e-05};
  double dTA11YSin[15]   = {0.0260681318,0.0287069673,0,0.02122356,0,0.00405945595,0.0130730001,0,0.00324590507,0,0.00718144537,0,0.00196385203,0,-5.79088536e-05};
  double dTA11YimSin[15] = {0,0,-0.0735281909,0,-0.0478521949,-0.00459199244,0,-0.0261676531,-0.00391062067,-5.5474199e-05,0,-0.0134519833,-0.00257691824,-7.54365475e-05,-5.79088536e-05};
  double dTA11ZSin[15]   = {0.292839946,0.210683976,0.0202988913,0.117572104,0.0173419978,0,0.0605086817,0.0113276397,0,0.00027136901,0.0301000943,0.00642577168,0,0.000327716963,5.33416943e-05};
  double dTA11ZimSin[15] = {0,0,0.0202988913,0,0.0173419978,0,0,0.0113276397,0,-0.00027136901,0,0.00642577168,0,-0.000327716963,0};
  
  double dTA22XSin[15]   = {0.0260681318,-0.0287069673,-0.0735281909,0.02122356,0.0478521949,-0.00405945595,-0.0130730001,-0.0261676531,0.00324590507,5.5474199e-05,0.00718144537,0.0134519833,-0.00196385203,-7.54365475e-05,-5.79088536e-05};
  double dTA22XimSin[15] = {0,0,0,0,0,-0.00459199244,0,0,0.00391062067,0,0,0,-0.00257691824,0,5.79088536e-05};
  double dTA22YSin[15]   = {0.0260681318,-0.0287069673,0,0.02122356,0,0.00405945595,-0.0130730001,0,-0.00324590507,0,0.00718144537,0,0.00196385203,0,-5.79088536e-05};
  double dTA22YimSin[15] = {0,0,-0.0735281909,0,0.0478521949,-0.00459199244,0,-0.0261676531,0.00391062067,-5.5474199e-05,0,0.0134519833,-0.00257691824,7.54365475e-05,-5.79088536e-05};
  double dTA22ZSin[15]   = {-0.292839946,0.210683976,-0.0202988913,-0.117572104,0.0173419978,0,0.0605086817,-0.0113276397,0,-0.00027136901,-0.0301000943,0.00642577168,0,0.000327716963,-5.33416943e-05};
  double dTA22ZimSin[15] = {0,0,-0.0202988913,0,0.0173419978,0,0,-0.0113276397,0,0.00027136901,0,0.00642577168,0,-0.000327716963,0};
  
  double dTA02XSin[15]   = {-0.0231772105,0.026943709,0.0650549213,-0.0205213805,-0.0448320512,0.00391889606,0.0128410777,0.0253052532,-0.00318379525,-5.53666116e-05,-0.00711276676,-0.0132282742,0.00194274847,7.53632455e-05,5.77965509e-05};
  double dTA02XimSin[15] = {0,0,0,0,0,0.00444733418,0,0,-0.00384593792,0,0,0,0.00255462452,0,-5.77965509e-05};
  double dTA02YSin[15]   = {-0.0231772105,0.026943709,0,-0.0205213805,0,-0.00391889606,0.0128410777,0,0.00318379525,0,-0.00711276676,0,-0.00194274847,0,5.77965509e-05};
  double dTA02YimSin[15] = {0,0,0.0650549213,0,-0.0448320512,0.00444733418,0,0.0253052532,-0.00384593792,5.53666116e-05,0,-0.0132282742,0.00255462452,-7.53632455e-05,5.77965509e-05};
  double dTA02ZSin[15]   = {0.233403613,-0.185886546,0.0190520793,0.110101501,-0.0167672253,0,-0.0585080455,0.011126415,0,0.000270318615,0.0295985635,-0.00636426116,0,-0.000327080644,5.32898619e-05};
  double dTA02ZimSin[15] = {0,0,0.0190520793,0,-0.0167672253,0,0,0.011126415,0,-0.000270318615,0,-0.00636426116,0,0.000327080644,0};
  
  double dTA10XSin[15]   = {-0.0229447351,0.026776716,0.0703936584,-0.0204764698,-0.049397592,0.0039097261,0.0128707059,0.0282749134,-0.00319178102,-5.53560674e-05,-0.00716797844,-0.0149071327,0.00196002174,7.53488931e-05,5.78319697e-05};
  double dTA10XimSin[15] = {0,0,0,0,0,0.00443828566,0,0,-0.00385420258,0,0,0,0.00257227267,0,-5.78319697e-05};
  double dTA10YSin[15]   = {-0.0229447351,0.026776716,0,-0.0204764698,0,-0.0039097261,0.0128707059,0,0.00319178102,0,-0.00716797844,0,-0.00196002174,0,5.78319697e-05};
  double dTA10YimSin[15] = {0,0,0.0703936584,0,-0.049397592,0.00443828566,0,0.0282749134,-0.00385420258,5.53560674e-05,0,-0.0149071327,0.00257227267,-7.53488931e-05,5.78319697e-05};
  double dTA10ZSin[15]   = {0.247041386,-0.201053506,0.0189339974,0.121301449,-0.0167304396,0,-0.0653717461,0.0111520383,0,0.00027043247,0.0333545636,-0.00641363872,0,-0.000327281034,5.32797133e-05};
  double dTA10ZimSin[15] = {0,0,0.0189339974,0,-0.0167304396,0,0,0.0111520383,0,-0.00027043247,0,-0.00641363872,0,0.000327281034,0};
  
  double dTA21XSin[15]   = {-0.00289092134,-0.00176325833,0.00847326962,-0.000702179466,0.00302014369,0.000140559883,-0.0002319224,0.000862399834,6.21098194e-05,-1.07587475e-07,-6.86786125e-05,0.000223709147,2.11035566e-05,-7.33020788e-08,1.12302729e-07};
  double dTA21XimSin[15] = {0,0,0,0,0,0.000144658265,0,0,6.46827515e-05,0,0,0,2.22937233e-05,0,-1.12302729e-07};
  double dTA21YSin[15]   = {-0.00289092134,-0.00176325833,0,-0.000702179466,0,-0.000140559883,-0.0002319224,0,-6.21098194e-05,0,-6.86786125e-05,0,-2.11035566e-05,0,1.12302729e-07};
  double dTA21YimSin[15] = {0,0,0.00847326962,0,0.00302014369,0.000144658265,0,0.000862399834,6.46827515e-05,1.07587475e-07,0,0.000223709147,2.22937233e-05,7.33020788e-08,1.12302729e-07};
  double dTA21ZSin[15]   = {-0.0594363331,-0.0247974305,-0.00124681193,-0.00747060322,-0.000574772552,0,-0.00200063617,-0.000201224723,0,-1.05039516e-06,-0.000501530811,-6.15105133e-05,0,-6.36319266e-07,-5.1832397e-08};
  double dTA21ZimSin[15] = {0,0,-0.00124681193,0,-0.000574772552,0,0,-0.000201224723,0,1.05039516e-06,0,-6.15105133e-05,0,6.36319266e-07,0};
  
  
  double dA00X[15]   = {2.06215379e-07,1.15146964e-05,-0.000224301532,8.95670998e-07,-2.52884084e-05,1.90441842e-08,-1.82688134e-08,-1.21983726e-06,4.50480896e-08,-2.73180646e-07,-1.26411065e-08,1.28380108e-08,9.87928444e-09,-4.93666544e-08,-1.63899461e-09};
  double dA00Y[15]   = {5.16775034e-06,0.000382210229,2.11820649e-06,5.64571239e-05,3.11430836e-07,1.47788249e-05,5.47993373e-06,1.32841834e-08,2.24921543e-06,2.09938378e-08,4.20475578e-07,-2.16817581e-09,2.38338795e-07,2.76421415e-09,1.83783143e-08};
  double dA00Z[15]   = {1.87271165e-05,0.000914832126,8.14694804e-06,8.27593478e-05,7.83466557e-07,3.81338114e-06,4.43664338e-06,-1.21374784e-08,6.40712427e-07,3.17674017e-08,1.02416474e-07,-1.11729335e-08,7.31267491e-08,6.27106198e-09,-2.14418796e-08};
  double dA00Xim[15] = {0,0,2.11695117e-06,0,3.11200437e-07,-6.97738829e-06,0,1.32590501e-08,-6.75884949e-07,1.5110907e-08,0,-2.17032212e-09,-2.37889085e-08,4.45754919e-09,-2.14455646e-08};
  double dA00Yim[15] = {0,0,-0.000258132337,0,-3.29153918e-05,-3.33372334e-07,0,-2.22534349e-06,-3.2762344e-08,2.13855916e-08,0,-8.9697776e-08,-1.82392461e-09,1.12796833e-08,-5.40467682e-10};
  double dA00Zim[15] = {0,0,0.000270268437,0,4.73901793e-05,3.11430836e-07,0,4.86660613e-06,1.74029773e-08,6.42311251e-07,0,3.84159121e-07,-3.01501753e-09,1.15022569e-07,2.55343314e-09};
  
  double dA11X[15]   = {-3.73051613e-05,0.000233890331,0.00195751593,0.000264963321,1.42474392e-05,-0.000358968663,3.56319725e-06,-4.08553733e-05,-1.09206199e-06,5.12263573e-05,-5.05655969e-06,9.14400233e-09,5.12097829e-06,2.25027579e-08,-6.56677879e-06};
  double dA11Y[15]   = {-1.59610808e-06,0.000179526803,-2.7008332e-05,-3.38278661e-05,-9.27260178e-06,1.11064453e-05,2.75654671e-06,2.31001838e-06,1.10066719e-06,-3.12209652e-06,-1.25482263e-07,-1.33503253e-07,-3.0339709e-07,-5.63701059e-08,5.49389424e-07};
  double dA11Z[15]   = {-9.75371433e-06,-0.000834208906,0.000165382798,-4.70809761e-05,0.000206760636,-2.5338756e-06,3.69756478e-05,3.15972105e-06,-3.06537793e-05,-1.14006184e-07,-1.0571562e-07,-4.48161513e-06,-4.83911517e-08,3.72586653e-06,7.00209075e-09};
  double dA11Xim[15] = {0,0,-2.71028716e-05,0,-9.27256367e-06,2.45392321e-05,0,2.31151329e-06,-3.73257583e-07,-4.33060915e-06,0,-1.33503076e-07,-2.26594772e-07,4.29290746e-08,6.53180024e-07};
  double dA11Yim[15] = {0,0,-0.00107574129,0,1.93153148e-05,0.000241216833,0,8.33158988e-06,-8.12813743e-07,-3.77837615e-05,0,7.69937578e-08,-1.9035557e-06,2.69937741e-09,5.17796552e-06};
  double dA11Zim[15] = {0,0,0.000126942288,0,-2.80270843e-05,-9.27245881e-06,0,2.44395824e-06,2.91592788e-06,2.96960949e-07,0,-1.16910372e-07,-1.90314466e-07,-2.81800636e-07,-4.75218193e-09};
  
  double dA22X[15]   = {3.39942321e-05,-0.00024542763,0.00168246196,-0.000210238617,1.13964319e-05,0.000293061468,-3.54429679e-06,-3.17637216e-05,1.04644438e-06,3.98952898e-05,3.82608209e-06,-2.17812233e-08,-3.82739369e-06,2.78138234e-08,4.86697278e-06};
  double dA22Y[15]   = {-1.51735318e-06,0.000170654521,2.48859655e-05,-3.18609107e-05,9.58419794e-06,1.11485522e-05,2.4798818e-06,-2.32325137e-06,1.08962595e-06,3.10094968e-06,-9.50887991e-08,1.31330149e-07,-3.01496718e-07,5.91363895e-08,5.47489465e-07};
  double dA22Z[15]   = {-9.45081085e-06,-0.000667329611,-0.000173545746,-4.25211224e-05,-0.000162856935,-3.59885902e-06,2.89784437e-05,-3.14704316e-06,-2.36347228e-05,8.19163544e-08,-6.4727528e-08,3.38785281e-06,-6.69692273e-08,-2.76484889e-06,7.34455968e-09};
  double dA22Xim[15] = {0,0,2.4981765e-05,0,9.58392985e-06,2.37979702e-05,0,-2.32472127e-06,-2.34530507e-07,4.31531708e-06,0,1.31327827e-07,-2.43785605e-07,-3.84656976e-08,6.54105716e-07};
  double dA22Yim[15] = {0,0,-0.000935241527,0,1.8593971e-05,-0.000198430565,0,6.24045835e-06,8.45795151e-07,-2.94275744e-05,0,7.20951619e-08,1.39140133e-06,7.04182202e-09,-3.82830186e-06};
  double dA22Zim[15] = {0,0,0.000120668638,0,-2.64015134e-05,9.58405497e-06,0,2.20158696e-06,-2.93326738e-06,3.49094126e-07,0,-8.94506447e-08,1.87292505e-07,-2.90004981e-07,7.30840512e-09};
  
  double dA01X[15]   = {-3.22345336e-06,-0.000236386116,0.000105435279,-3.44887049e-05,1.12234674e-05,7.62763375e-06,-3.29207548e-06,4.38354836e-07,1.0765013e-06,1.5530967e-07,-2.46687853e-07,-2.19076227e-08,1.03709381e-07,2.7375113e-08,-2.34342145e-09};
  double dA01Y[15]   = {-2.54437284e-06,-0.000186662351,-4.59241476e-05,-2.72445528e-05,-1.02940521e-05,-7.36760497e-06,-2.60161144e-06,-1.34724049e-06,-1.11902864e-06,-3.75311258e-07,-1.95041557e-07,-1.36130971e-07,-1.18217179e-07,-5.83224088e-08,-1.01265141e-08};
  double dA01Z[15]   = {-9.21162624e-06,-0.000444151179,-0.000166672703,-3.91008227e-05,-2.88854015e-05,-2.43755564e-06,-1.97019789e-06,-2.9175438e-06,-4.06165378e-07,1.61711328e-07,-3.06924622e-08,-2.24926542e-07,-4.5873934e-08,2.87252405e-08,1.09100635e-08};
  double dA01Xim[15] = {0,0,-4.60529004e-05,0,-1.03223602e-05,3.11930258e-06,0,-1.35094269e-06,2.68692865e-07,2.99969261e-07,0,-1.36510225e-07,3.30397728e-09,4.29035917e-08,1.11995815e-08};
  double dA01Yim[15] = {0,0,0.000127179928,0,1.60959576e-05,5.62953427e-06,0,1.07555843e-06,6.77409079e-07,4.25257641e-09,0,4.24051613e-08,4.98342276e-08,-3.48792385e-09,-7.60563582e-09};
  double dA01Zim[15] = {0,0,-0.000131992711,0,-2.28818388e-05,-1.02933868e-05,0,-2.31209545e-06,-1.71708871e-06,-3.47179082e-07,0,-1.7835007e-07,-1.94201191e-07,-6.16110309e-08,-5.4911821e-09};
  
  double dA12X[15]   = {3.42135593e-05,-6.27321838e-07,-0.00207832845,-0.000300454086,1.24165839e-07,0.000366545745,2.01881353e-08,4.00620575e-05,-1.883732e-08,-5.13633655e-05,4.76718766e-06,-2.88725047e-09,-5.00927387e-06,2.56470654e-09,6.56440714e-06};
  double dA12Y[15]   = {-7.82133474e-07,8.66617433e-08,6.71649796e-05,6.51284704e-06,-1.37172113e-08,-1.81532429e-05,-1.76468226e-09,-1.04946481e-06,1.50004478e-09,3.4712356e-06,-9.72901843e-08,1.94391778e-10,1.82993527e-07,-1.67283531e-10,-5.58740241e-07};
  double dA12Z[15]   = {-1.06513985e-09,0.0013256407,-1.11393185e-07,-7.77410958e-08,-0.000236567742,3.38656688e-08,-3.41310332e-05,1.32703167e-08,3.10855601e-05,-5.47161647e-09,2.08244839e-09,4.2184094e-06,-1.81378718e-09,-3.70361199e-06,7.34397488e-10};
  double dA12Xim[15] = {0,0,6.71765149e-05,0,2.36686894e-09,-2.10153946e-05,0,-1.05044989e-06,-1.98301653e-10,4.04541043e-06,0,1.44899053e-11,2.45094168e-07,2.1734907e-11,-6.43378847e-07};
  double dA12Yim[15] = {0,0,0.000933435666,0,-2.79847531e-08,-0.000235073519,0,-9.79283414e-06,6.71499629e-09,3.77526108e-05,0,4.17885627e-10,1.97244781e-06,-1.14800944e-09,-5.18347226e-06};
  double dA12Zim[15] = {0,0,5.93703362e-08,0,5.12382743e-06,-1.35422553e-08,0,-2.44689453e-10,-1.30909098e-06,1.9523326e-09,0,-8.62015805e-08,5.99993126e-11,2.27121203e-07,-2.52521719e-10};
  
  double dA20X[15]   = {-3.29817109e-06,0.000244789903,0.000107326171,-3.63850367e-05,-1.15231754e-05,7.55582105e-06,3.56492725e-06,4.49928864e-07,-1.06568675e-06,1.74149213e-07,-2.76741212e-07,2.47316907e-08,1.01842255e-07,-3.04353745e-08,-7.67522001e-10};
  double dA20Y[15]   = {2.29939773e-06,-0.000170743294,4.2274197e-05,2.53484743e-05,-9.59845833e-06,7.00365491e-06,-2.47802636e-06,1.27385874e-06,-1.09120375e-06,3.70102977e-07,1.9237056e-07,-1.3112518e-07,1.18515279e-07,-5.93124311e-08,1.12226885e-08};
  double dA20Z[15]   = {9.45153398e-06,-0.000464936929,0.000173415825,4.26021938e-05,-3.05924848e-05,3.56172411e-06,-2.34788044e-06,3.16090102e-06,-6.03834723e-07,-8.78515513e-08,6.2552294e-08,-2.51983907e-07,6.88653322e-08,1.5949582e-08,-8.13825246e-09};
  double dA20Xim[15] = {0,0,4.2192837e-05,0,-9.58153778e-06,-2.78324579e-06,0,1.27430392e-06,2.34740676e-07,-2.70035668e-07,0,-1.31313803e-07,-1.30148115e-09,3.84887256e-08,-1.07478222e-08};
  double dA20Yim[15] = {0,0,0.000138547392,0,-1.85653164e-05,6.47731531e-06,0,1.38708726e-06,-8.38886385e-07,6.11280519e-08,0,-7.25270054e-08,7.07144173e-08,-5.85927417e-09,-4.9580178e-09};
  double dA20Zim[15] = {0,0,-0.000120729501,0,2.12776989e-05,-9.59814041e-06,0,-2.20129522e-06,1.62421877e-06,-3.51138385e-07,0,1.75646498e-07,-1.87226467e-07,6.28928575e-08,-7.57331769e-09};
  
  double dA300SingX[15]  = {9.86198779e-05,0,-0.12410435,0.0450956638,0,-0.00944337877,0,-0.0721386116,0,0.000151368039,0.0185932465,0,-0.00559566867,0,-0.000168690941};
  double dA300SingY[15]  = {9.86198779e-05,0,0,0.0450956638,0,0.00944337877,0,0,0,0,0.0185932465,0,0.00559566867,0,-0.000168690941};
  double dA300SingZ[15]  = {0,0.411918429,0,0,0.0438523272,0,0.205282586,0,0,0,0,0.0202705215,0,0.000957234103,0};
  double dA300SingImX[15] = {0,0,0,0,0,-0.0107127037,0,0,0,0,0,0,-0.00735635572,0,0.000168690941};
  double dA300SingImY[15] = {0,0,-0.12410435,0,0,-0.0107127037,0,-0.0721386116,0,-0.000151368039,0,0,-0.00735635572,0,-0.000168690941};
  double dA300SingImZ[15] = {0,0,0,0,0.0438523272,0,0,0,0,0,0,0.0202705215,0,-0.000957234103,0};
  double dA301SingX[15]  = {-4.9309939e-05,0.0220829716,0.062052175,-0.0225478319,-0.0562133905,0.00472168939,0.0157698682,0.0360693058,-0.00429877247,-7.56840197e-05,-0.00929662327,-0.0196499761,0.00279783433,0.000109650465,8.43454704e-05};
  double dA301SingY[15]  = {-4.9309939e-05,0.0220829716,0,-0.0225478319,0,-0.00472168939,0.0157698682,0,0.00429877247,0,-0.00929662327,0,-0.00279783433,0,8.43454704e-05};
  double dA301SingZ[15]  = {0.000669040352,-0.205959215,0.018909177,0.169025954,-0.0219261636,0,-0.102641293,0.0165399661,0,0.000371781162,0.0574897253,-0.0101352607,0,-0.000478617052,7.77358922e-05};
  double dA301SingImX[15] = {0,0,0,0,0,0.00535635185,0,0,-0.00519207948,0,0,0,0.00367817786,0,-8.43454704e-05};
  double dA301SingImY[15] = {0,0,0.062052175,0,-0.0562133905,0.00535635185,0,0.0360693058,-0.00519207948,7.56840197e-05,0,-0.0196499761,0.00367817786,-0.000109650465,8.43454704e-05};
  double dA301SingImZ[15] = {0,0,0.018909177,0,-0.0219261636,0,0,0.0165399661,0,-0.000371781162,0,-0.0101352607,0,0.000478617052,0};
  double dA311SingX[15]  = {5.72648594e-05,0.0240810393,-0.0645550478,0.023556662,-0.0546181652,-0.00490462892,0.0161710219,-0.0328472471,-0.00437673347,7.58461701e-05,0.00935395855,-0.0176235661,-0.00280306604,0.00010977891,-8.44576108e-05};
  double dA311SingY[15]  = {5.72648594e-05,0.0240810393,0,0.023556662,0,0.00490462892,0.0161710219,0,0.00437673347,0,0.00935395855,0,0.00280306604,0,-8.44576108e-05};
  double dA311SingZ[15]  = {0.000777188964,0.216006242,0.0200614413,0.163450533,0.0226948821,0,0.0967182298,0.0166934257,0,0.000373060839,0.052503108,0.0101614797,0,0.000479253036,7.78262717e-05};
  double dA311SingImX[15] = {0,0,0,0,0,-0.00554411518,0,0,-0.00527335242,0,0,0,-0.00368455502,0,8.44576108e-05};
  double dA311SingImY[15] = {0,0,-0.0645550478,0,-0.0546181652,-0.00554411518,0,-0.0328472471,-0.00527335242,-7.58461701e-05,0,-0.0176235661,-0.00368455502,-0.00010977891,-8.44576108e-05};
  double dA311SingImZ[15] = {0,0,0.0200614413,0,0.0226948821,0,0,0.0166934257,0,-0.000373060839,0,0.0101614797,0,-0.000479253036,0};
  double dA312SingX[15]  = {-2.80117153e-06,-0.000676765177,0.00850510006,-0.000310770918,0.00422992022,0.000186929527,-8.9293014e-05,0.00144395057,0.000100269823,-1.01035007e-07,-2.53367642e-05,0.000416253951,4.23779999e-05,-3.42732322e-08,1.79889046e-07};
  double dA312SingY[15]  = {-2.80117153e-06,-0.000676765177,0,-0.000310770918,0,-0.000186929527,-8.9293014e-05,0,-0.000100269823,0,-2.53367642e-05,0,-4.23779999e-05,0,1.79889046e-07};
  double dA312SingZ[15]  = {-0.000129018856,-0.0174514139,-0.00149893734,-0.00501352187,-0.000904353254,0,-0.000599983973,-0.000340929889,0,-1.29429374e-06,0.000175912376,-9.33831912e-05,0,-6.69695205e-07,-9.04773458e-08};
  double dA312SingImX[15] = {0,0,0,0,0,0.000194098223,0,0,0.000106424551,0,0,0,4.653186e-05,0,-1.79889046e-07};
  double dA312SingImY[15] = {0,0,0.00850510006,0,0.00422992022,0.000194098223,0,0.00144395057,0.000106424551,1.01035007e-07,0,0.000416253951,4.653186e-05,3.42732322e-08,1.79889046e-07};
  double dA312SingImZ[15] = {0,0,-0.00149893734,0,-0.000904353254,0,0,-0.000340929889,0,1.29429374e-06,0,-9.33831912e-05,0,6.69695205e-07,0};
  double dA322SingX[15]  = {5.72648594e-05,-0.0240810393,-0.0645550478,0.023556662,0.0546181652,-0.00490462892,-0.0161710219,-0.0328472471,0.00437673347,7.58461701e-05,0.00935395855,0.0176235661,-0.00280306604,-0.00010977891,-8.44576108e-05};
  double dA322SingY[15]  = {5.72648594e-05,-0.0240810393,0,0.023556662,0,0.00490462892,-0.0161710219,0,-0.00437673347,0,0.00935395855,0,0.00280306604,0,-8.44576108e-05};
  double dA322SingZ[15]  = {-0.000777188964,0.216006242,-0.0200614413,-0.163450533,0.0226948821,0,0.0967182298,-0.0166934257,0,-0.000373060839,-0.052503108,0.0101614797,0,0.000479253036,-7.78262717e-05};
  double dA322SingImX[15] = {0,0,0,0,0,-0.00554411518,0,0,0.00527335242,0,0,0,-0.00368455502,0,8.44576108e-05};
  double dA322SingImY[15] = {0,0,-0.0645550478,0,0.0546181652,-0.00554411518,0,-0.0328472471,0.00527335242,-7.58461701e-05,0,0.0176235661,-0.00368455502,0.00010977891,-8.44576108e-05};
  double dA322SingImZ[15] = {0,0,-0.0200614413,0,0.0226948821,0,0,-0.0166934257,0,0.000373060839,0,0.0101614797,0,-0.000479253036,0};
  double dA320SingX[15]  = {-5.44636878e-05,0.0234042741,0.0560499477,-0.0232458911,-0.050388245,0.0047176994,0.0160817289,0.0314032966,-0.00427646365,-7.57451351e-05,-0.00932862179,-0.0172073122,0.00276068804,0.000109744637,8.42777217e-05};
  double dA320SingY[15]  = {-5.44636878e-05,0.0234042741,0,-0.0232458911,0,-0.0047176994,0.0160817289,0,0.00427646365,0,-0.00932862179,0,-0.00276068804,0,8.42777217e-05};
  double dA320SingZ[15]  = {0.000648170108,-0.198554828,0.018562504,0.158437011,-0.0217905289,0,-0.0961182459,0.0163524958,0,0.000371766545,0.0526790204,-0.0100680966,0,-0.000478583341,7.77357943e-05};
  double dA320SingImX[15] = {0,0,0,0,0,0.00535001695,0,0,-0.00516692786,0,0,0,0.00363802317,0,-8.42777217e-05};
  double dA320SingImY[15] = {0,0,0.0560499477,0,-0.050388245,0.00535001695,0,0.0314032966,-0.00516692786,7.57451351e-05,0,-0.0172073122,0.00363802317,-0.000109744637,8.42777217e-05};
  double dA320SingImZ[15] = {0,0,0.018562504,0,-0.0217905289,0,0,0.0163524958,0,-0.000371766545,0,-0.0100680966,0,0.000478583341,0};

  double L0[15] = {0.0676056845,0.00750146809,0.00185347754,0.000437801593,
    0.000208024957,1.36224131e-05,2.00421651e-05,1.64224844e-05,1.78797695e-06,
    -5.89301269e-07,6.82090982e-07,1.08352568e-06,1.70746019e-07,
    -8.90524479e-08,-2.94499438e-08};
  double L0_im[15] = {0,0,0.0014651663,0,0.000164472954,5.73054024e-05,0,
    1.29876352e-05,7.52793049e-06,1.2692253e-06,0,8.57198339e-07,
    7.19654983e-07,1.91648298e-07,1.49038432e-08};
  double L1[15] = {0.0698331075,-0.00794092114,-0.0017799131,0.000487357134,
    0.000205679316,1.33067114e-05,-2.43860909e-05,-1.69189375e-05,
    -1.80158025e-06,4.36309276e-07,1.00017433e-06,1.17946738e-06,
    1.78774053e-07,-6.80727362e-08,-2.14820881e-08};
  double L1_im[15] = {0,0,-0.00136887958,0,0.00015821673,4.99704248e-05,0,
    -1.30187435e-05,-6.77152753e-06,-1.04925268e-06,0,9.07947072e-07,
    6.72700267e-07,1.63546868e-07,1.23857389e-08};
  
  double L0Sing[15] = {2.12314203,0,0.0648975111,0.400806897,0,0,0,0.0236506803,0,0.000386005325,0.0653401901,0,0,0,5.91781253e-05};
  double L1Sing[15] = {1.47781839,0.585679893,0.0368659055,0.209779027,0.0234391407,0,0.0782857646,0.0122571145,0,0.000194454835,0.0302378247,0.00584777108,0,0.000205135689,2.96522067e-05};
  double L2Sing[15] = {1.47781839,-0.585679893,0.0368659055,0.209779027,-0.0234391407,0,-0.0782857646,0.0122571145,0,0.000194454835,0.0302378247,-0.00584777108,0,-0.000205135689,2.96522067e-05};
  double L0SingIm[15] = {0,0,0.0648975111,0,0,0,0,0.0236506803,0,-0.000386005325,0,0,0,0,0};
  double L1SingIm[15] = {0,0,0.0368659055,0,0.0234391407,0,0,0.0122571145,0,-0.000194454835,0,0.00584777108,0,-0.000205135689,0};
  double L2SingIm[15] = {0,0,0.0368659055,0,-0.0234391407,0,0,0.0122571145,0,-0.000194454835,0,-0.00584777108,0,0.000205135689,0};
  
  double dLdx0[15]   = {8.40101049e-05,1.24331395e-05,-0.000242192425,7.31717594e-07,-2.06593419e-05,1.55581287e-08,-1.33241252e-08,-8.89672694e-07,3.28552475e-08,-1.99240807e-07,-8.66579646e-09,8.80077931e-09,6.77249795e-09,-3.38420832e-08,-1.12357203e-09};
  double dLdx0im[15] = {0,0,2.28580488e-06,0,2.5423491e-07,-5.70017092e-06,0,9.67031849e-09,-4.92948036e-07,1.10209466e-08,0,-1.4878104e-09,-1.63078951e-08,3.05576208e-09,-1.47014739e-08};
  double dLdy0[15]   = {0.00210529035,0.000412696344,2.28716034e-06,4.61225952e-05,2.54423134e-07,1.20735473e-05,3.99671952e-06,9.68864912e-09,1.64043648e-06,1.53115868e-08,2.88246584e-07,-1.48633905e-09,1.63387238e-07,1.89493832e-09,1.25987967e-08};
  double dLdy0im[15] = {0,0,-0.000278721666,0,-2.68901989e-05,-2.72348221e-07,0,-1.62302579e-06,-2.38947962e-08,1.55973075e-08,0,-6.14900814e-08,-1.25034619e-09,7.73250663e-09,-3.70504189e-10};
  double dLdz0[15]   = {0.00762924195,0.000987801596,8.79677052e-06,6.76101727e-05,6.4005228e-07,3.1153382e-06,3.23580906e-06,-8.85231453e-09,4.67295407e-07,2.31691478e-08,7.02090686e-08,-7.65932691e-09,5.01302257e-08,4.29897071e-09,-1.46989478e-08};
  double dLdz0im[15] = {0,0,0.000291825774,0,3.87153632e-05,2.54423134e-07,0,3.54939688e-06,1.26926387e-08,4.68461488e-07,0,2.6335074e-07,-2.06687034e-09,7.8850864e-08,1.75044264e-09};
  
  double dLdx1[15]   = {-0.0151977536,0.000252546053,0.00211365266,0.000216461541,1.16394323e-05,-0.000293259117,2.5987723e-06,-2.97973435e-05,-7.96481432e-07,3.73612879e-05,-3.46639885e-06,6.26844357e-09,3.51055942e-06,1.54262065e-08,-4.50169202e-06};
  double dLdx1im[15] = {0,0,-2.92646695e-05,0,-7.57521232e-06,2.0047303e-05,0,1.68587262e-06,-2.7223064e-07,-3.15847435e-06,0,-9.1519716e-08,-1.55336415e-07,2.94289603e-08,4.47771334e-07};
  double dLdy1[15]   = {-0.000650238638,0.000193846344,-2.91625891e-05,-2.7635644e-05,-7.57524345e-06,9.07340022e-06,2.01045206e-06,1.68478234e-06,8.02757528e-07,-2.27706113e-06,-8.60212474e-08,-9.15198376e-08,-2.07986336e-07,-3.86431253e-08,3.76620268e-07};
  double dLdy1im[15] = {0,0,-0.0011615453,0,1.57796286e-05,0.000197061869,0,6.07653842e-06,-5.92815299e-07,-2.75571028e-05,0,5.27811574e-08,-1.30493531e-06,1.8504911e-09,3.54962559e-06};
  double dLdz1[15]   = {-0.0039735667,-0.000900747651,0.000178574175,-3.84627599e-05,0.00016891291,-2.07004733e-06,2.69677154e-05,2.3045021e-06,-2.23569415e-05,-8.31489506e-08,-7.2470716e-08,-3.07225989e-06,-3.31733516e-08,2.55417522e-06,4.80010931e-09};
  double dLdz1im[15] = {0,0,0.000137067547,0,-2.28967006e-05,-7.57512666e-06,0,1.78246966e-06,2.12669467e-06,2.16584667e-07,0,-8.01450002e-08,-1.30465353e-07,-1.93181424e-07,-3.25774024e-09};
  
  double dLdx2[15]   = {0.0138489138,-0.000265003597,0.00181665965,-0.000171754244,9.31030448e-06,0.000239416296,-2.58498749e-06,-2.31664638e-05,7.63210809e-07,2.90971189e-05,2.6228755e-06,-1.49315764e-08,-2.62377464e-06,1.90670755e-08,3.33643225e-06};
  double dLdx2im[15] = {0,0,2.69743778e-05,0,7.82958263e-06,1.94417298e-05,0,-1.6955057e-06,-1.71051823e-07,3.14732127e-06,0,9.00285282e-08,-1.67121164e-07,-2.63692031e-08,4.4840592e-07};
  double dLdy2[15]   = {-0.000618154669,0.000184266385,2.6870937e-05,-2.60287416e-05,7.82980165e-06,9.10779938e-06,1.80867005e-06,-1.69443365e-06,7.94704739e-07,2.26163795e-06,-6.51857634e-08,9.00301203e-08,-2.06683583e-07,4.05394823e-08,3.75317798e-07};
  double dLdy2im[15] = {0,0,-0.00100983891,0,1.51903274e-05,-0.000162107666,0,4.55139841e-06,6.16869867e-07,-2.14626247e-05,0,4.94230467e-08,9.53840604e-07,4.8273461e-09,-2.62439721e-06};
  double dLdz2[15]   = {-0.00385016682,-0.000720557614,-0.000187388221,-3.47375917e-05,-0.000133045823,-2.94008455e-06,2.11350569e-05,-2.29525565e-06,-1.72376826e-05,5.97446444e-08,-4.43723484e-08,2.32245831e-06,-4.59090896e-08,-1.89537345e-06,5.03488039e-09};
  double dLdz2im[15] = {0,0,0.000130293493,0,-2.15686919e-05,7.82968485e-06,0,1.60569927e-06,-2.13934102e-06,2.54607332e-07,0,-6.13206667e-08,1.28393723e-07,-1.98805708e-07,5.01009553e-09};
  
  double dLdx0Sing[15]   = {0.0401767624,0,-0.134003247,0.0368408609,0,-0.00771475956,0,-0.0526133729,0,0.000110398342,0.0127461381,0,-0.00383597161,0,-0.000115641883};
  double dLdy0Sing[15]   = {0.0401767624,0,0,0.0368408609,0,0.00771475956,0,0,0,0,0.0127461381,0,0.00383597161,0,-0.000115641883};
  double dLdz0Sing[15]   = {0,0.444774151,0,0,0.0358251182,0,0.149720226,0,0,0,0,0.0138959523,0,0.000656208053,0};
  double dLdx0imSing[15] = {0,0,0,0,0,-0.00875173338,0,0,0,0,0,0,-0.00504296686,0,0.000115641883};
  double dLdy0imSing[15] = {0,0,-0.134003247,0,0,-0.00875173338,0,-0.0526133729,0,-0.000110398342,0,0,-0.00504296686,0,-0.000115641883};
  double dLdz0imSing[15] = {0,0,0,0,0.0358251182,0,0,0,0,0,0,0.0138959523,0,-0.000656208053,0};
  double dLdx1Sing[15]   = {0.023329137,0.0260018077,-0.069704132,0.0192445933,-0.0446202596,-0.00400683207,0.0117941278,-0.0239567192,-0.00319211453,5.53174332e-05,0.00641237385,-0.0120813978,-0.00192157227,7.52562041e-05,-5.78978164e-05};
  double dLdy1Sing[15]   = {0.023329137,0.0260018077,0,0.0192445933,0,0.00400683207,0.0117941278,0,0.00319211453,0,0.00641237385,0,0.00192157227,0,-5.78978164e-05};
  double dLdz1Sing[15]   = {0.316619093,0.233235481,0.0216615958,0.133530762,0.018540563,0,0.0705402028,0.0121751364,0,0.000272087147,0.0359922011,0.00696594994,0,0.00032854001,5.33518666e-05};
  double dLdx1imSing[15] = {0,0,0,0,0,-0.00452925977,0,0,-0.00384605208,0,0,0,-0.00252585514,0,5.78978164e-05};
  double dLdy1imSing[15] = {0,0,-0.069704132,0,-0.0446202596,-0.00452925977,0,-0.0239567192,-0.00384605208,-5.53174332e-05,0,-0.0120813978,-0.00252585514,-7.52562041e-05,-5.78978164e-05};
  double dLdz1imSing[15] = {0,0,0.0216615958,0,0.018540563,0,0,0.0121751364,0,-0.000272087147,0,0.00696594994,0,-0.00032854001,0};

} ; // end ASolverUTest


TEST_F(ASolverUTest, checkGamma)
{
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());

  EXPECT_NEAR( ASolvTest.get_gamma_ni( 0, 1).real(),  1.463995711, preclim);
  EXPECT_NEAR( ASolvTest.get_gamma_ni( 0, 5).real(),  1.760111936, preclim);
  
  EXPECT_NEAR( ASolvTest.get_gamma_ni( 1, 2).real(),  1.621243794, preclim);
  EXPECT_NEAR( ASolvTest.get_gamma_ni( 1, 7).real(),  1.799701878, preclim);
}

TEST_F(ASolverUTest, checkDelta)
{
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  
  EXPECT_NEAR( ASolvTest.get_delta_ni( 0, 1).real(),   0.87554313, preclim);
  EXPECT_NEAR( ASolvTest.get_delta_ni( 0, 5).real(),   0.06832297, preclim);
  
  EXPECT_NEAR( ASolvTest.get_delta_ni( 1, 2).real(),   11.4370663, preclim);
  EXPECT_NEAR( ASolvTest.get_delta_ni( 1, 7).real()/181.9847, 1.0, preclim);
}


TEST_F(ASolverUTest, checkE)
{
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  
  EXPECT_NEAR(ASolvTest.get_E_ni( 0, 0, 0).real(), 5.0, preclim);
  EXPECT_NEAR(ASolvTest.get_E_ni( 0, 0, 0).imag(), 0.0, preclim);

  EXPECT_NEAR(ASolvTest.get_E_ni( 0, 5, 0).real()/-4.7683716e-06, 1.0, preclim);
  EXPECT_NEAR(ASolvTest.get_E_ni( 0, 5, 0).imag(),                0.0, preclim);

  EXPECT_NEAR(ASolvTest.get_E_ni( 0, 6, -5).real(),        0.0, preclim);
  EXPECT_NEAR(ASolvTest.get_E_ni( 0, 6, -5).imag(),        0.0, preclim);

  EXPECT_NEAR(ASolvTest.get_E_ni( 1, 0, 0).real(),-0.4, preclim);
  EXPECT_NEAR(ASolvTest.get_E_ni( 1, 0, 0).imag(), 0.0, preclim);

  EXPECT_NEAR(ASolvTest.get_E_ni( 1, 3, -2).real(), 0.00014227949, preclim);
  EXPECT_NEAR(ASolvTest.get_E_ni( 1, 3, -2).imag(), 0.00134791098, preclim);

  EXPECT_NEAR(ASolvTest.get_E_ni( 1, 6, -5).real()/-6.7122985e-06, 1., preclim);
  EXPECT_NEAR(ASolvTest.get_E_ni( 1, 6, -5).imag()/-3.8653666e-06, 1., preclim);
}

TEST_F(ASolverUTest, checkSH)
{
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());

  EXPECT_NEAR(ASolvTest.get_SH_ij( 0, 0, 0, 0).real(), 1.0, preclim);
  EXPECT_NEAR(ASolvTest.get_SH_ij( 0, 0, 0, 0).imag(), 0.0, preclim);

  EXPECT_NEAR(ASolvTest.get_SH_ij( 0, 0, 5, 0).real(),-1.0, preclim);
  EXPECT_NEAR(ASolvTest.get_SH_ij( 0, 0, 5, 0).imag(), 0.0, preclim);

  EXPECT_NEAR(ASolvTest.get_SH_ij( 0, 0, 6,-5).real(), 0.0, preclim);
  EXPECT_NEAR(ASolvTest.get_SH_ij( 0, 0, 6,-5).imag(), 0.0, preclim);
  
  EXPECT_NEAR(ASolvTest.get_SH_ij( 1, 0, 0, 0).real(), 1.0, preclim);
  EXPECT_NEAR(ASolvTest.get_SH_ij( 1, 0, 0, 0).imag(), 0.0, preclim);

  EXPECT_NEAR(ASolvTest.get_SH_ij( 1, 0, 3,-2).real(),-0.0522108211, preclim);
  EXPECT_NEAR(ASolvTest.get_SH_ij( 1, 0, 3,-2).imag(), 0.4946288312, preclim);
  
  EXPECT_NEAR(ASolvTest.get_SH_ij( 1, 0, 6, 5).real(), 0.3615494298, preclim);
  EXPECT_NEAR(ASolvTest.get_SH_ij( 1, 0, 6, 5).imag(),-0.2082030597, preclim);
}

TEST_F(ASolverUTest, checkAMulti)
{
  mol_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = { Pt(0.0,0.0,-5.0), Pt(10.0,7.8,25.0), Pt(-10.0,7.8,25.0)};
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                    molInd, 0);
    mol_.push_back( molNew );
  }
  
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000);
  
  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      EXPECT_NEAR(ASolvTest.get_A_ni(0, n, m).real()/ATrip0[ct], 1.0, preclim);
      EXPECT_NEAR(ASolvTest.get_A_ni(1, n, m).real()/ATrip1[ct], 1.0, preclim);
      EXPECT_NEAR(ASolvTest.get_A_ni(2, n, m).real()/ATrip2[ct], 1.0, preclim);
      
      if (ATrip0im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(0, n, m).imag()/ATrip0im[ct],
                    1.0, preclim);
      if (ATrip1im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(1, n, m).imag()/ATrip1im[ct],
                    1.0, preclim);
      if (ATrip2im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(2, n, m).imag()/ATrip2im[ct],
                    1.0, preclim);
      ct++;
    }
  }
}


TEST_F(ASolverUTest, checkAMultiPBC)
{
  mol_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[4] = { Pt(0.0,0.0,-5.0), Pt(10.0,7.8,25.0),
                Pt(-10.0,7.8,25.0), Pt( 100, 100, 100)};
  for (int molInd = 0; molInd < 4; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_.push_back( molNew );
  }
  
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_, 85, 175);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000);
  
  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      EXPECT_NEAR(ASolvTest.get_A_ni(0, n, m).real()/ATrip0[ct], 1.0, preclim);
      EXPECT_NEAR(ASolvTest.get_A_ni(1, n, m).real()/ATrip1[ct], 1.0, preclim);
      EXPECT_NEAR(ASolvTest.get_A_ni(2, n, m).real()/ATrip2[ct], 1.0, preclim);
      
      if (ATripPBC3[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(3, n, m).real()/ATripPBC3[ct],
                    1.0,preclim);
      
      if (ATrip0im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(0, n, m).imag()/ATrip0im[ct],
                    1.0, preclim);
      if (ATrip1im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(1, n, m).imag()/ATrip1im[ct],
                    1.0, preclim);
      if (ATrip2im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(2, n, m).imag()/ATrip2im[ct],
                    1.0, preclim);
      if (ATripPBC3im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(3, n, m).imag()/ATripPBC3im[ct],
                    1.0, preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkAPBC)
{
  mol_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = { Pt(0.0,0.0,0.0), Pt(10.0,7.8,25.0), Pt(-10.0,-7.8,-25.0)};
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_.push_back( molNew );
  }
  
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_, 40, 175);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000);

  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      EXPECT_NEAR(1.0, 1.0, preclim);
      EXPECT_NEAR(ASolvTest.get_A_ni(0,n,m).real()/APBCTrip0[ct],1.0,preclim);
      EXPECT_NEAR(ASolvTest.get_A_ni(1,n,m).real()/APBCTrip1[ct],1.0,preclim);
      EXPECT_NEAR(ASolvTest.get_A_ni(2,n,m).real()/APBCTrip2[ct],1.0,preclim);
      
      if (APBCTrip0im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(0, n, m).imag()/APBCTrip0im[ct],
                    1.0, preclim);
      if (APBCTrip1im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(1, n, m).imag()/APBCTrip1im[ct],
                    1.0, preclim);
      if (APBCTrip2im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(2, n, m).imag()/APBCTrip2im[ct],
                    1.0, preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkAPBCoutside)
{
  mol_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = { Pt(0.0,0.0,0.0), Pt(10.0,7.8,25.0), Pt(-10.0,-7.8,-25.0)};
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_.push_back( molNew );
  }
  
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_, 25, 150);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000);
  
  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      if (APBCoutside[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(0, n, m).real()/APBCoutside[ct],
                    1.0, preclim);
      if (APBCoutside[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(1, n, m).real()/APBCoutside[ct],
                    1.0, preclim);
      if (APBCoutside[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(2, n, m).real()/APBCoutside[ct],
                    1.0, preclim);

      if (APBCoutsideim[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(0, n, m).imag()/APBCoutsideim[ct],
                    1.0, preclim);
      if (APBCoutsideim[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(1, n, m).imag()/APBCoutsideim[ct],
                    1.0, preclim);
      if (APBCoutsideim[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(2, n, m).imag()/APBCoutsideim[ct],
                    1.0, preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkA)
{
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-20, 1000);

  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      EXPECT_NEAR(ASolvTest.get_A_ni( 0, n, m).real()/A0[ct],      1, preclim);
      if (A0_im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni( 0, n, m).imag()/A0_im[ct], 1, preclim);
      EXPECT_NEAR(ASolvTest.get_A_ni( 1, n, m).real()/A1[ct],      1, preclim);
      if (A1_im[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni( 1, n, m).imag()/A1_im[ct], 1, preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkASingMult)
{
  mol_sing_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[2] = {  Pt( 0.0, 0.0, -5.0 ), Pt( 0.0, 0.0, 0.0 )};
  for (int molInd = 0; molInd < 2; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_sing_.push_back( molNew );
  }
  
  const int vals           = 5;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_sing_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000);
  
  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      if (AMul0Sing[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(0,n,m).real()/AMul0Sing[ct], 1, preclim);
      if (AMul0Singim[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(0,n,m).imag()/AMul0Singim[ct],
                    1, preclim);
      if (AMul1Sing[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(1,n,m).real()/AMul1Sing[ct], 1, preclim);
      if (AMul1Singim[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(1,n,m).imag()/AMul1Singim[ct],
                    1, preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkASingMultFlip)
{
  mol_sing_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = {  Pt( 0.0, 0.0, 0.0 ),Pt( 0.0, 0.0, -5.0 ),Pt( 0.0, 0.0, 5.0)};
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_sing_.push_back( molNew );
  }
  
  const int vals           = 5;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_sing_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000);
  
  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      if (AMul0SingF[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(0,n,m).real()/AMul0SingF[ct], 1, preclim);
      if (AMul0SingimF[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(0,n,m).imag()/AMul0SingimF[ct],
                    1, preclim);
      if (AMul1SingF[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(1,n,m).real()/AMul1SingF[ct], 1, preclim);
      if (AMul1SingimF[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(1,n,m).imag()/AMul1SingimF[ct],
                    1, preclim);
      
      if (AMul2SingF[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(2,n,m).real()/AMul2SingF[ct], 1, preclim);
      if (AMul2SingimF[ct] != 0)
        EXPECT_NEAR(ASolvTest.get_A_ni(2,n,m).imag()/AMul2SingimF[ct],
                    1, preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkASing)
{
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_sing_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-30, 1000);
  
  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      EXPECT_NEAR( ASolvTest.get_A_ni( 0, n, m).real(), A0Sing[ct], preclim);
      EXPECT_NEAR( ASolvTest.get_A_ni( 0, n, m).imag(),        0.0, preclim);
      EXPECT_NEAR( ASolvTest.get_A_ni( 1, n, m).real(), A1Sing[ct], preclim);
      EXPECT_NEAR( ASolvTest.get_A_ni( 1, n, m).imag(),        0.0, preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkgradT_A)
{
  mol_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = { Pt( 0.0, 0.0, -5.0 ),
    Pt( 10.0, 7.8, 25.0 ),Pt(-10.0, 7.8, 25.0) };
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_.push_back( molNew );
  }
  
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000); ASolvTest.solve_gradA(1E-30, 1000);
  
  VecOfMats<cmplx>::type dT_A00 = ASolvTest.get_gradT_Aij( 0, 0);
  VecOfMats<cmplx>::type dT_A11 = ASolvTest.get_gradT_Aij( 1, 1);
  VecOfMats<cmplx>::type dT_A22 = ASolvTest.get_gradT_Aij( 2, 2);
  
  VecOfMats<cmplx>::type dT_A02 = ASolvTest.get_gradT_Aij( 0, 2);
  VecOfMats<cmplx>::type dT_A10 = ASolvTest.get_gradT_Aij( 1, 0);
  VecOfMats<cmplx>::type dT_A21 = ASolvTest.get_gradT_Aij( 2, 1);
  
  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      if (dTA00X[ct] != 0) EXPECT_NEAR(dT_A00[0](n, m+vals).real()/dTA00X[ct],
                                       1.0, preclim);
      if (dTA00Y[ct] != 0) EXPECT_NEAR(dT_A00[1](n, m+vals).real()/dTA00Y[ct],
                                       1.0, preclim);
      if (dTA00Z[ct] != 0) EXPECT_NEAR(dT_A00[2](n, m+vals).real()/dTA00Z[ct],
                                       1.0, preclim);
      if (dTA00Xim[ct] != 0) EXPECT_NEAR(dT_A00[0](n, m+vals).imag()/
                                         dTA00Xim[ct], 1.0, preclim);
      if (dTA00Yim[ct] != 0) EXPECT_NEAR(dT_A00[1](n, m+vals).imag()/
                                         dTA00Yim[ct], 1.0, preclim);
      if (dTA00Zim[ct] != 0) EXPECT_NEAR(dT_A00[2](n, m+vals).imag()/
                                         dTA00Zim[ct], 1.0, preclim);
      
      if (dTA11X[ct] != 0) EXPECT_NEAR(dT_A11[0](n, m+vals).real()/dTA11X[ct],
                                       1.0, preclim);
      if (dTA11Y[ct] != 0) EXPECT_NEAR(dT_A11[1](n, m+vals).real()/dTA11Y[ct],
                                       1.0, preclim);
      if (dTA11Z[ct] != 0) EXPECT_NEAR(dT_A11[2](n, m+vals).real()/dTA11Z[ct],
                                       1.0, preclim);
      if (dTA11Xim[ct] != 0) EXPECT_NEAR(dT_A11[0](n, m+vals).imag()/
                                         dTA11Xim[ct], 1.0, preclim);
      if (dTA11Yim[ct] != 0) EXPECT_NEAR(dT_A11[1](n, m+vals).imag()/
                                         dTA11Yim[ct], 1.0, preclim);
      if (dTA11Zim[ct] != 0) EXPECT_NEAR(dT_A11[2](n, m+vals).imag()/
                                         dTA11Zim[ct], 1.0, preclim);
      
      if (dTA22X[ct] != 0) EXPECT_NEAR(dT_A22[0](n, m+vals).real()/dTA22X[ct],
                                       1.0, preclim);
      if (dTA22Y[ct] != 0) EXPECT_NEAR(dT_A22[1](n, m+vals).real()/dTA22Y[ct],
                                       1.0, preclim);
      if (dTA22Z[ct] != 0) EXPECT_NEAR(dT_A22[2](n, m+vals).real()/dTA22Z[ct],
                                       1.0, preclim);
      if (dTA22Xim[ct] != 0) EXPECT_NEAR(dT_A22[0](n, m+vals).imag()/
                                         dTA22Xim[ct], 1.0, preclim);
      if (dTA22Yim[ct] != 0) EXPECT_NEAR(dT_A22[1](n, m+vals).imag()/
                                         dTA22Yim[ct], 1.0, preclim);
      if (dTA22Zim[ct] != 0) EXPECT_NEAR(dT_A22[2](n, m+vals).imag()/
                                         dTA22Zim[ct], 1.0, preclim);
      
      if (dTA02X[ct] != 0) EXPECT_NEAR(dT_A02[0](n, m+vals).real()/dTA02X[ct],
                                       1.0, preclim);
      if (dTA02Y[ct] != 0) EXPECT_NEAR(dT_A02[1](n, m+vals).real()/dTA02Y[ct],
                                       1.0, preclim);
      if (dTA02Z[ct] != 0) EXPECT_NEAR(dT_A02[2](n, m+vals).real()/dTA02Z[ct],
                                       1.0, preclim);
      if (dTA02Xim[ct] != 0) EXPECT_NEAR(dT_A02[0](n, m+vals).imag()/
                                         dTA02Xim[ct], 1.0, preclim);
      if (dTA02Yim[ct] != 0) EXPECT_NEAR(dT_A02[1](n, m+vals).imag()/
                                         dTA02Yim[ct], 1.0, preclim);
      if (dTA02Zim[ct] != 0) EXPECT_NEAR(dT_A02[2](n, m+vals).imag()/
                                         dTA02Zim[ct], 1.0, preclim);
      
      if (dTA10X[ct] != 0) EXPECT_NEAR(dT_A10[0](n, m+vals).real()/dTA10X[ct],
                                       1.0, preclim);
      if (dTA10Y[ct] != 0) EXPECT_NEAR(dT_A10[1](n, m+vals).real()/dTA10Y[ct],
                                       1.0, preclim);
      if (dTA10Z[ct] != 0) EXPECT_NEAR(dT_A10[2](n, m+vals).real()/dTA10Z[ct],
                                       1.0, preclim);
      if (dTA10Xim[ct] != 0) EXPECT_NEAR(dT_A10[0](n, m+vals).imag()/
                                         dTA10Xim[ct], 1.0, preclim);
      if (dTA10Yim[ct] != 0) EXPECT_NEAR(dT_A10[1](n, m+vals).imag()/
                                         dTA10Yim[ct], 1.0, preclim);
      if (dTA10Zim[ct] != 0) EXPECT_NEAR(dT_A10[2](n, m+vals).imag()/
                                         dTA10Zim[ct], 1.0, preclim);
      
      if (dTA21X[ct] != 0) EXPECT_NEAR(dT_A21[0](n, m+vals).real()/dTA21X[ct],
                                       1.0, preclim);
      if (dTA21Y[ct] != 0) EXPECT_NEAR(dT_A21[1](n, m+vals).real()/dTA21Y[ct],
                                       1.0, preclim);
      if (dTA21Z[ct] != 0) EXPECT_NEAR(dT_A21[2](n, m+vals).real()/dTA21Z[ct],
                                       1.0, preclim);
      if (dTA21Xim[ct] != 0) EXPECT_NEAR(dT_A21[0](n, m+vals).imag()/
                                         dTA21Xim[ct], 1.0, preclim);
      if (dTA21Yim[ct] != 0) EXPECT_NEAR(dT_A21[1](n, m+vals).imag()/
                                         dTA21Yim[ct], 1.0, preclim);
      if (dTA21Zim[ct] != 0) EXPECT_NEAR(dT_A21[2](n, m+vals).imag()/
                                         dTA21Zim[ct], 1.0, preclim);
      ct++;
    }
  }
}


TEST_F(ASolverUTest, checkdT_ASingFlip)
{
  mol_sing_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = {  Pt( 0.0, 0.0, 0.0 ),Pt( 0.0, 0.0, -5.0 ),Pt( 0.0, 0.0, 5.0)};
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_sing_.push_back( molNew );
  }
  const int vals = 5;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_sing_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000); ASolvTest.solve_gradA(1E-30, 1000);
  
  VecOfMats<cmplx>::type dT_A00 = ASolvTest.get_gradT_Aij( 0, 0);
  VecOfMats<cmplx>::type dT_A11 = ASolvTest.get_gradT_Aij( 1, 1);
  VecOfMats<cmplx>::type dT_A22 = ASolvTest.get_gradT_Aij( 2, 2);
  
  VecOfMats<cmplx>::type dT_A02 = ASolvTest.get_gradT_Aij( 0, 2);
  VecOfMats<cmplx>::type dT_A10 = ASolvTest.get_gradT_Aij( 1, 0);
  VecOfMats<cmplx>::type dT_A21 = ASolvTest.get_gradT_Aij( 2, 1);

  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      if (dTA00XSin[ct] != 0) EXPECT_NEAR(dT_A00[0](n, m+vals).real()/
                                          dTA00XSin[ct], 1.0, preclim);
      if (dTA00YSin[ct] != 0) EXPECT_NEAR(dT_A00[1](n, m+vals).real()/
                                          dTA00YSin[ct], 1.0, preclim);
      if (dTA00ZSin[ct] != 0) EXPECT_NEAR(dT_A00[2](n, m+vals).real()/
                                          dTA00ZSin[ct], 1.0, preclim);
      if (dTA00XimSin[ct] != 0) EXPECT_NEAR(dT_A00[0](n, m+vals).imag()/
                                         dTA00XimSin[ct], 1.0, preclim);
      if (dTA00YimSin[ct] != 0) EXPECT_NEAR(dT_A00[1](n, m+vals).imag()/
                                         dTA00YimSin[ct], 1.0, preclim);
      if (dTA00ZimSin[ct] != 0) EXPECT_NEAR(dT_A00[2](n, m+vals).imag()/
                                         dTA00ZimSin[ct], 1.0, preclim);
      
      if (dTA11XSin[ct] != 0) EXPECT_NEAR(dT_A11[0](n, m+vals).real()/
                                          dTA11XSin[ct], 1.0, preclim);
      if (dTA11YSin[ct] != 0) EXPECT_NEAR(dT_A11[1](n, m+vals).real()/
                                          dTA11YSin[ct], 1.0, preclim);
      if (dTA11ZSin[ct] != 0) EXPECT_NEAR(dT_A11[2](n, m+vals).real()/
                                          dTA11ZSin[ct], 1.0, preclim);
      if (dTA11XimSin[ct] != 0) EXPECT_NEAR(dT_A11[0](n, m+vals).imag()/
                                         dTA11XimSin[ct], 1.0, preclim);
      if (dTA11YimSin[ct] != 0) EXPECT_NEAR(dT_A11[1](n, m+vals).imag()/
                                         dTA11YimSin[ct], 1.0, preclim);
      if (dTA11ZimSin[ct] != 0) EXPECT_NEAR(dT_A11[2](n, m+vals).imag()/
                                         dTA11ZimSin[ct], 1.0, preclim);
      
      if (dTA22XSin[ct] != 0) EXPECT_NEAR(dT_A22[0](n, m+vals).real()/
                                          dTA22XSin[ct], 1.0, preclim);
      if (dTA22YSin[ct] != 0) EXPECT_NEAR(dT_A22[1](n, m+vals).real()/
                                          dTA22YSin[ct], 1.0, preclim);
      if (dTA22ZSin[ct] != 0) EXPECT_NEAR(dT_A22[2](n, m+vals).real()/
                                          dTA22ZSin[ct], 1.0, preclim);
      if (dTA22XimSin[ct] != 0) EXPECT_NEAR(dT_A22[0](n, m+vals).imag()/
                                         dTA22XimSin[ct], 1.0, preclim);
      if (dTA22YimSin[ct] != 0) EXPECT_NEAR(dT_A22[1](n, m+vals).imag()/
                                         dTA22YimSin[ct], 1.0, preclim);
      if (dTA22ZimSin[ct] != 0) EXPECT_NEAR(dT_A22[2](n, m+vals).imag()/
                                         dTA22ZimSin[ct], 1.0, preclim);
      
      if (dTA02XSin[ct] != 0) EXPECT_NEAR(dT_A02[0](n, m+vals).real()/
                                          dTA02XSin[ct], 1.0, preclim);
      if (dTA02YSin[ct] != 0) EXPECT_NEAR(dT_A02[1](n, m+vals).real()/
                                          dTA02YSin[ct], 1.0, preclim);
      if (dTA02ZSin[ct] != 0) EXPECT_NEAR(dT_A02[2](n, m+vals).real()/
                                          dTA02ZSin[ct], 1.0, preclim);
      if (dTA02XimSin[ct] != 0) EXPECT_NEAR(dT_A02[0](n, m+vals).imag()/
                                         dTA02XimSin[ct], 1.0, preclim);
      if (dTA02YimSin[ct] != 0) EXPECT_NEAR(dT_A02[1](n, m+vals).imag()/
                                         dTA02YimSin[ct], 1.0, preclim);
      if (dTA02ZimSin[ct] != 0) EXPECT_NEAR(dT_A02[2](n, m+vals).imag()/
                                         dTA02ZimSin[ct], 1.0, preclim);
      
      if (dTA10XSin[ct] != 0) EXPECT_NEAR(dT_A10[0](n, m+vals).real()/
                                          dTA10XSin[ct], 1.0, preclim);
      if (dTA10YSin[ct] != 0) EXPECT_NEAR(dT_A10[1](n, m+vals).real()/
                                          dTA10YSin[ct], 1.0, preclim);
      if (dTA10ZSin[ct] != 0) EXPECT_NEAR(dT_A10[2](n, m+vals).real()/
                                          dTA10ZSin[ct], 1.0, preclim);
      if (dTA10XimSin[ct] != 0) EXPECT_NEAR(dT_A10[0](n, m+vals).imag()/
                                         dTA10XimSin[ct], 1.0, preclim);
      if (dTA10YimSin[ct] != 0) EXPECT_NEAR(dT_A10[1](n, m+vals).imag()/
                                         dTA10YimSin[ct], 1.0, preclim);
      if (dTA10ZimSin[ct] != 0) EXPECT_NEAR(dT_A10[2](n, m+vals).imag()/
                                         dTA10ZimSin[ct], 1.0, preclim);
      
      if (dTA21XSin[ct] != 0) EXPECT_NEAR(dT_A21[0](n, m+vals).real()/
                                          dTA21XSin[ct], 1.0, preclim);
      if (dTA21YSin[ct] != 0) EXPECT_NEAR(dT_A21[1](n, m+vals).real()/
                                          dTA21YSin[ct], 1.0, preclim);
      if (dTA21ZSin[ct] != 0) EXPECT_NEAR(dT_A21[2](n, m+vals).real()/
                                          dTA21ZSin[ct], 1.0, preclim);
      if (dTA21XimSin[ct] != 0) EXPECT_NEAR(dT_A21[0](n, m+vals).imag()/
                                         dTA21XimSin[ct], 1.0, preclim);
      if (dTA21YimSin[ct] != 0) EXPECT_NEAR(dT_A21[1](n, m+vals).imag()/
                                         dTA21YimSin[ct], 1.0, preclim);
      if (dTA21ZimSin[ct] != 0) EXPECT_NEAR(dT_A21[2](n, m+vals).imag()/
                                         dTA21ZimSin[ct], 1.0, preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkgradA)
{
  mol_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = { Pt(0.0,0.0,-5.0), Pt(10.0,7.8,25.0), Pt(-10.0,7.8,25.0) };
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_.push_back( molNew );
  }
  
  const int vals = 5;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-25, 1000); ASolvTest.solve_gradA(1E-50, 1000);
  
  int ct = 0;
  for ( int n = 0; n < 3; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      if (dA00X[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(0,0,n,m).real()
                                       /dA00X[ct], 1.0, preclim);
      if (dA00Y[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(0,0,n,m).real()
                                      /dA00Y[ct], 1.0, preclim);
      if (dA00Z[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(0,0,n,m).real()
                                      /dA00Z[ct], 1.0, preclim);
      if (dA00Xim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(0,0,n,m).imag()
                                        /dA00Xim[ct], 1.0, preclim);
      if (dA00Yim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(0,0,n,m).imag()
                                        /dA00Yim[ct], 1.0, preclim);
      if (dA00Zim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(0,0,n,m).imag()
                                        /dA00Zim[ct], 1.0, preclim);
      
      if (dA11X[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(1,1,n,m).real()
                                      /dA11X[ct], 1.0, preclim);
      if (dA11Y[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(1,1,n,m).real()
                                      /dA11Y[ct], 1.0, preclim);
      if (dA11Z[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(1,1,n,m).real()
                                      /dA11Z[ct], 1.0, preclim);
      if (dA11Xim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(1,1,n,m).imag()
                                        /dA11Xim[ct], 1.0, preclim);
      if (dA11Yim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(1,1,n,m).imag()
                                        /dA11Yim[ct], 1.0, preclim);
      if (dA11Zim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(1,1,n,m).imag()
                                        /dA11Zim[ct], 1.0, preclim);
      
      if (dA22X[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(2,2,n,m).real()
                                      /dA22X[ct], 1.0, preclim);
      if (dA22Y[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(2,2,n,m).real()
                                      /dA22Y[ct], 1.0, preclim);
      if (dA22Z[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(2,2,n,m).real()
                                      /dA22Z[ct], 1.0, preclim);
      if (dA22Xim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(2,2,n,m).imag()
                                        /dA22Xim[ct], 1.0, preclim);
      if (dA22Yim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(2,2,n,m).imag()
                                        /dA22Yim[ct], 1.0, preclim);
      if (dA22Zim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(2,2,n,m).imag()
                                        /dA22Zim[ct], 1.0, preclim);
      
      if (dA01X[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(0,1,n,m).real()/
                                      dA01X[ct], 1.0, preclim);
      if (dA01Y[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(0,1,n,m).real()/
                                      dA01Y[ct], 1.0, preclim);
      if (dA01Z[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(0,1,n,m).real()/
                                      dA01Z[ct], 1.0, preclim);
      if (dA01Xim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(0,1,n,m).imag()/
                                         dA01Xim[ct], 1.0, preclim);
      if (dA01Yim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(0,1,n,m).imag()/
                                         dA01Yim[ct], 1.0, preclim);
      if (dA01Zim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(0,1,n,m).imag()/
                                         dA01Zim[ct], 1.0, preclim);
      
      if (dA12X[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(1,2,n,m).real()/
                                      dA12X[ct], 1.0, preclim);
      if (dA12Y[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(1,2,n,m).real()/
                                      dA12Y[ct], 1.0, preclim);
      if (dA12Z[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(1,2,n,m).real()/
                                      dA12Z[ct], 1.0, preclim);
      if (dA12Xim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(1,2,n,m).imag()/
                                        dA12Xim[ct], 1.0, preclim);
      if (dA12Yim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(1,2,n,m).imag()/
                                        dA12Yim[ct], 1.0, preclim);
      if (dA12Zim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(1,2,n,m).imag()/
                                        dA12Zim[ct], 1.0, preclim);
      
      if (dA20X[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(2,0,n,m).real()/
                                      dA20X[ct], 1.0, preclim);
      if (dA20Y[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(2,0,n,m).real()/
                                      dA20Y[ct], 1.0, preclim);
      if (dA20Z[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(2,0,n,m).real()/
                                      dA20Z[ct], 1.0, preclim);
      if (dA20Xim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdx_ni(2,0,n,m).imag()/
                                        dA20Xim[ct], 1.0, preclim);
      if (dA20Yim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdy_ni(2,0,n,m).imag()/
                                        dA20Yim[ct], 1.0, preclim);
      if (dA20Zim[ct] != 0) EXPECT_NEAR(ASolvTest.get_dAdz_ni(2,0,n,m).imag()/
                                        dA20Zim[ct], 1.0, preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkgradASing)
{
  mol_sing_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = {  Pt( 0.0, 0.0, 0.0 ),Pt( 0.0, 0.0, -5.0 ),Pt( 0.0, 0.0, 5.0)};
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_sing_.push_back( molNew );
  }
  const int vals           = 5;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_sing_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000); ASolvTest.solve_gradA(1E-30, 1000);

  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      if ( dA300SingX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 0, 0, n, m).real()/dA300SingX[ct],
                  1, preclim);
      if ( dA300SingY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 0, 0, n, m).real()/dA300SingY[ct],
                    1, preclim);
      if ( dA300SingZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 0, 0, n, m).real()/dA300SingZ[ct],
                    1, preclim);
      if ( dA300SingImX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 0, 0, n, m).imag()/dA300SingImX[ct],
                    1, preclim);
      if ( dA300SingImY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 0, 0, n, m).imag()/dA300SingImY[ct],
                    1, preclim);
      if ( dA300SingImZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 0, 0, n, m).imag()/dA300SingImZ[ct],
                    1, preclim);
      
      if ( dA301SingX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 0, 1, n, m).real()/dA301SingX[ct],
                    1, preclim);
      if ( dA301SingY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 0, 1, n, m).real()/dA301SingY[ct],
                    1, preclim);
      if ( dA301SingZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 0, 1, n, m).real()/dA301SingZ[ct],
                    1, preclim);
      if ( dA301SingImX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 0, 1, n, m).imag()/dA301SingImX[ct],
                    1, preclim);
      if ( dA301SingImY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 0, 1, n, m).imag()/dA301SingImY[ct],
                    1, preclim);
      if ( dA301SingImZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 0, 1, n, m).imag()/dA301SingImZ[ct],
                    1, preclim);

      if ( dA311SingX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 1, 1, n, m).real()/dA311SingX[ct],
                    1, preclim);
      if ( dA311SingY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 1, 1, n, m).real()/dA311SingY[ct],
                    1, preclim);
      if ( dA311SingZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 1, 1, n, m).real()/dA311SingZ[ct],
                    1, preclim);
      if ( dA311SingImX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 1, 1, n, m).imag()/dA311SingImX[ct],
                    1, preclim);
      if ( dA311SingImY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 1, 1, n, m).imag()/dA311SingImY[ct],
                    1, preclim);
      if ( dA311SingImZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 1, 1, n, m).imag()/dA311SingImZ[ct],
                    1, preclim);
      
      if ( dA312SingX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 1, 2, n, m).real()/dA312SingX[ct],
                    1, preclim);
      if ( dA312SingY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 1, 2, n, m).real()/dA312SingY[ct],
                    1, preclim);
      if ( dA312SingZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 1, 2, n, m).real()/dA312SingZ[ct],
                    1, preclim);
      if ( dA312SingImX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 1, 2, n, m).imag()/dA312SingImX[ct],
                    1, preclim);
      if ( dA312SingImY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 1, 2, n, m).imag()/dA312SingImY[ct],
                    1, preclim);
      if ( dA312SingImZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 1, 2, n, m).imag()/dA312SingImZ[ct],
                    1, preclim);
      
      if ( dA322SingX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 2, 2, n, m).real()/dA322SingX[ct],
                    1, preclim);
      if ( dA322SingY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 2, 2, n, m).real()/dA322SingY[ct],
                    1, preclim);
      if ( dA322SingZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 2, 2, n, m).real()/dA322SingZ[ct],
                    1, preclim);
      if ( dA322SingImX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 2, 2, n, m).imag()/dA322SingImX[ct],
                    1, preclim);
      if ( dA322SingImY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 2, 2, n, m).imag()/dA322SingImY[ct],
                    1, preclim);
      if ( dA322SingImZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 2, 2, n, m).imag()/dA322SingImZ[ct],
                    1, preclim);
      
      if ( dA320SingX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 2, 0, n, m).real()/dA320SingX[ct],
                    1, preclim);
      if ( dA320SingY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 2, 0, n, m).real()/dA320SingY[ct],
                    1, preclim);
      if ( dA320SingZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 2, 0, n, m).real()/dA320SingZ[ct],
                    1, preclim);
      if ( dA322SingImX[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdx_ni( 2, 2, n, m).imag()/dA322SingImX[ct],
                    1, preclim);
      if ( dA322SingImY[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdy_ni( 2, 2, n, m).imag()/dA322SingImY[ct],
                    1, preclim);
      if ( dA322SingImZ[ct] != 0 )
        EXPECT_NEAR( ASolvTest.get_dAdz_ni( 2, 2, n, m).imag()/dA322SingImZ[ct],
                    1, preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkL)
{
  mol_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = { Pt( 0.0, 0.0, -5.0 ),
    Pt( 10.0, 7.8, 25.0 ),Pt(-10.0, 7.8, 25.0) };
  for (int molInd = 0; molInd < 2; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_.push_back( molNew );
  }
  const int vals = nvals;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000);
  shared_ptr<VecOfMats<cmplx>::type> myL = ASolvTest.get_L();
  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      EXPECT_NEAR( myL->operator[](0)(n,m+vals).real()/L0[ct], 1.0, preclim);
      if (L0_im[ct] != 0)
        EXPECT_NEAR(myL->operator[](0)(n,m+vals).imag()/L0_im[ct],1.0,preclim);
      EXPECT_NEAR( myL->operator[](1)(n,m+vals).real()/L1[ct], 1.0, preclim);
      if (L1_im[ct] != 0)
        EXPECT_NEAR(myL->operator[](1)(n,m+vals).imag()/L1_im[ct],1.0,preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkLSing)
{
  mol_sing_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = {  Pt( 0.0, 0.0, 0.0 ),Pt( 0.0, 0.0, -5.0 ),Pt( 0.0, 0.0, 5.0)};
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_sing_.push_back( molNew );
  }
  const int vals = 5;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_sing_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000);
  shared_ptr<VecOfMats<cmplx>::type> myL = ASolvTest.get_L();
  
  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      if (L0Sing[ct] != 0) EXPECT_NEAR( myL->operator[](0)(n,m+vals).real()/
                                       L0Sing[ct], 1,preclim);
      if (L0SingIm[ct] != 0)
        EXPECT_NEAR(myL->operator[](0)(n,m+vals).imag()/L0SingIm[ct],1,preclim);
      if (L1Sing[ct] != 0)
        EXPECT_NEAR(myL->operator[](1)(n,m+vals).real()/L1Sing[ct], 1,preclim);
      if (L1SingIm[ct] != 0)
        EXPECT_NEAR(myL->operator[](1)(n,m+vals).imag()/L1SingIm[ct],
                    1, preclim);
      if (L2Sing[ct] != 0)
        EXPECT_NEAR(myL->operator[](2)(n,m+vals).real()/L2Sing[ct], 1,preclim);
      if (L2SingIm[ct] != 0)
        EXPECT_NEAR(myL->operator[](2)(n,m+vals).imag()/L2SingIm[ct],1,preclim);
      ct++;
    }
  }
}

TEST_F(ASolverUTest, checkdL)
{
  mol_.clear( );
  shared_ptr<MoleculeAM> molNew;
  Pt pos[3] = { Pt( 0.0, 0.0, -5.0 ),
  Pt( 10.0, 7.8, 25.0 ),Pt(-10.0, 7.8, 25.0) };
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0] = pos[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1] = pos[molInd] + Pt(1.0, 0.0, 0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2] = pos[molInd] + Pt(0.0, 1.0, 0.0);
    
    molNew = make_shared<MoleculeAM> ( "stat", 2.0, charges, posCharges, vdW, pos[molInd],
                                      molInd, 0);
    mol_.push_back( molNew );
  }
  
  const int vals           = 5;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-20, 1000); ASolvTest.solve_gradA(1E-30, 1000);
  shared_ptr<MyVector<VecOfMats<cmplx>::type > > mydL = ASolvTest.get_gradL();

  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      if (dLdx0[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[0](n,m+vals).real()/dLdx0[ct],
                    1, preclim);
      if (dLdx0im[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[0](n,m+vals).imag()/dLdx0im[ct],
                    1, preclim);
      if (dLdy0[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[1](n,m+vals).real()/dLdy0[ct],
                    1, preclim);
      if (dLdy0im[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[1](n,m+vals).imag()/dLdy0im[ct],
                    1, preclim);
      if (dLdz0[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[2](n,m+vals).real()/dLdz0[ct],
                    1, preclim);
      if (dLdz0im[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[2](n,m+vals).imag()/dLdz0im[ct],
                    1, preclim);
      
      if (dLdx1[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[0](n,m+vals).real()/dLdx1[ct],
                    1, preclim);
      if (dLdx1im[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[0](n,m+vals).imag()/dLdx1im[ct],
                    1, preclim);
      if (dLdy1[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[1](n,m+vals).real()/dLdy1[ct],
                    1, preclim);
      if (dLdy1im[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[1](n,m+vals).imag()/dLdy1im[ct],
                    1, preclim);
      if (dLdz1[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[2](n,m+vals).real()/dLdz1[ct],
                    1, preclim);
      if (dLdz1im[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[2](n,m+vals).imag()/dLdz1im[ct],
                    1, preclim);
      
      if (dLdx2[ct] != 0)
        EXPECT_NEAR( mydL->operator[](2)[0](n,m+vals).real()/dLdx2[ct],
                    1, preclim);
      if (dLdx2im[ct] != 0)
        EXPECT_NEAR( mydL->operator[](2)[0](n,m+vals).imag()/dLdx2im[ct],
                    1, preclim);
      if (dLdy2[ct] != 0)
        EXPECT_NEAR( mydL->operator[](2)[1](n,m+vals).real()/dLdy2[ct],
                    1, preclim);
      if (dLdy2im[ct] != 0)
        EXPECT_NEAR( mydL->operator[](2)[1](n,m+vals).imag()/dLdy2im[ct],
                    1, preclim);
      if (dLdz2[ct] != 0)
        EXPECT_NEAR( mydL->operator[](2)[2](n,m+vals).real()/dLdz2[ct],
                    1, preclim);
      if (dLdz2im[ct] != 0)
        EXPECT_NEAR( mydL->operator[](2)[2](n,m+vals).imag()/dLdz2im[ct],
                    1, preclim);
      ct++;
    }
  }  
}

TEST_F(ASolverUTest, checkdLSing)
{
  mol_sing_.clear();
  shared_ptr<MoleculeAM> molSing;
  Pt cgPosSi[3] = { Pt( 0.0, 0.0, 0.0 ), Pt( 0.0, 0.0, -5.0 ),
    Pt( 0.0, 0.0, 5.0 )};
  
  for (int molInd = 0; molInd < 3; molInd ++ )
  {
    int M = 3; vector<double> charges(M); vector<double> vdW(M);
    vector<Pt> posCharges(M);
    charges[0]=2.0; vdW[0]=0; posCharges[0]=cgPosSi[molInd];
    charges[1]=2.0; vdW[1]=0; posCharges[1]=cgPosSi[molInd] + Pt(1.0,0.0,0.0);
    charges[2]=2.0; vdW[2]=0; posCharges[2]=cgPosSi[molInd] + Pt(0.0,1.0,0.0);
    
    molSing = make_shared<MoleculeAM>( "stat", 2.0, charges, posCharges, vdW, cgPosSi[molInd],
                     molInd, 0);
    mol_sing_.push_back( molSing );
  }
  const int vals           = 5;
  shared_ptr<BesselConstants> bConsta = make_shared<BesselConstants>(2*vals);
  shared_ptr<BesselCalc> bCalcu = make_shared<BesselCalc>(2*vals, bConsta);
  shared_ptr<SHCalcConstants> SHConsta = make_shared<SHCalcConstants>(2*vals);
  shared_ptr<SHCalc> SHCalcu = make_shared<SHCalc>(2*vals, SHConsta);
  shared_ptr<SystemAM> sys = make_shared<SystemAM>(mol_sing_);
  
  ASolver ASolvTest(bCalcu, SHCalcu, sys, const_, vals, sys->get_cutoff());
  ASolvTest.solve_A(1E-40, 1000); ASolvTest.solve_gradA(1E-40, 1000);
  shared_ptr<MyVector<VecOfMats<cmplx>::type > > mydL = ASolvTest.get_gradL();
  
  int ct = 0;
  for ( int n = 0; n < 5; n++ )
  {
    for ( int m = 0; m <= n; m++ )
    {
      if (dLdx0Sing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[0](n,m+vals).real()/dLdx0Sing[ct],
                    1, preclim);
      if (dLdy0Sing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[1](n,m+vals).real()/dLdy0Sing[ct],
                    1, preclim);
      if (dLdz0Sing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[2](n,m+vals).real()/dLdz0Sing[ct],
                    1, preclim);
      if (dLdx0imSing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[0](n,m+vals).imag()/dLdx0imSing[ct],
                    1, preclim);
      if (dLdy0imSing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[1](n,m+vals).imag()/dLdy0imSing[ct],
                    1, preclim);
      if (dLdz0imSing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](0)[2](n,m+vals).imag()/dLdz0imSing[ct],
                    1, preclim);
      if (dLdx1Sing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[0](n,m+vals).real()/dLdx1Sing[ct],
                    1, preclim);
      if (dLdy1Sing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[1](n,m+vals).real()/dLdy1Sing[ct],
                    1, preclim);
      if (dLdz1Sing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[2](n,m+vals).real()/dLdz1Sing[ct],
                    1, preclim);
      if (dLdx1imSing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[0](n,m+vals).imag()/dLdx1imSing[ct],
                    1, preclim);
      if (dLdy1imSing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[1](n,m+vals).imag()/dLdy1imSing[ct],
                    1, preclim);
      if (dLdz1imSing[ct] != 0)
        EXPECT_NEAR( mydL->operator[](1)[2](n,m+vals).imag()/dLdz1imSing[ct],
                    1, preclim);
      ct++;
    }
  }
}


#endif
